diff --git a/private/Synthesizer/Basic/NumberTheory.hs b/private/Synthesizer/Basic/NumberTheory.hs new file mode 100644 --- /dev/null +++ b/private/Synthesizer/Basic/NumberTheory.hs @@ -0,0 +1,896 @@ +{-# LANGUAGE NoImplicitPrelude #-} +{- +Some of these functions might be moved to NumericPrelude. + +Wikipedia: (primitive) roots of unity modulo n + (primitive) roots must be units and all units are (primitive) roots + maximum possible order for primitive roots - Carmichael + all possible orders: divisor of Carmichael (proof? statement already in Carmichael-function-article) + sum of primitive roots that vanishes + order of primitive root is a divisor of each possible exponent + proof with GCD and diophantine in exponent + check for primitive root: fast exponentiation, + primitivity: check exponents that are prime divisors + how to find a primitive root: just try + sum of powers of a primitive root is zero + number of primitive roots of given order + g(n,k) > 0 if k|lambda(n) + g(n,k) = 0 else + g(n,1) = 1 + g(4,2) = 1 + g(2^n,2) = 3 for n>=3 ((-1) is always a square root of 1) + g(2^n,2^k) = 2^k for k>=2 && k=3 and n in OEIS:A033948 + sum(g(n,k), k\in\N) = phi(n) + There are only a few patterns that occur as rows of g, + but a row of g (i.e. g(n)) does functionally depend on + either lambda(n) nor phi(n) + lambda(14) = 6 nozeros(g(14)) = [1,1,2,2] (f ~ [1,2,3,6]) + lambda(21) = 6 nozeros(g(21)) = [1,3,2,6] (f ~ [1,4,3,12]) + phi(13) = 12 nozeros(g(13)) = [1,1,2,2,2,4] (f ~ [1,2,3,4,6,12]) + phi(21) = 12 nozeros(g(21)) = [1,3,2,6] (f ~ [1,4,3,12]) + However length(nozeros(f(n))) = numberofdivisors(lambda(n)) + numberofdivisors=A000005 + number of roots of given order + easier to compute + k|m => f(n,k) | f(n,m) + g(n,k) = f(n,k) - sum(f(n,d), d|k and k/d prime) + ... + inclusion-exclusion-principle + better to write the other round: + f(n,k) = sum(g(n,d), d|k) - this is Dirichlet convolution + RUNM says f(n,k) is multiplicative + list it in multiplicative function + f(n,1) = 1 for n>=2 + f(n,lambda(n)) = phi(n) + f(n,a·b) = f(n,a)·f(n,b) if a and b are coprime (conjecture) + f(n,lcm(a,b)) = lcm(f(n,a),f(n,b)) (conjecture) + If this conjecture is true, then we only need to know f(n,p^i). + The following conjecture is wrong: + for prime p it is f(n,p^i) = gcd(lambda(n),p^i) + counterexamples + f(8,2) = 4, lambda(8)=2 + f(63,3) = 9, lambda(63)=6 + f(275,5) = 25, lambda(275)=20 + f(1247,7) = 49, lambda(1247)=84 + It seems to be: + for prime p it is f(n,p^i) = p^j for some j + How to find a modulus where there is a primitive root of order o? + just try numbers from the sequence o+1, 2*o+1, 3*o+1 + Because of [[Dirichlet's theorem on arithmetic progressions]] + you will somewhen find a prime p, + and its Carmichael value is p-1, which is a multiple of o. + In this ring even 'o' is a unit. + How to find a modulus where there are primitive roots of orders o1,..,ok? + Just search for a modulus with roots of order lcm(o1,...,ok). + The smallest such modulus should also be the smallest one + that has primitive roots of orders o1,..,ok? + Proof: If a ring has primitive roots of orders o1,..,ok + then all orders divide the carmichael value of that ring, + thus lcm(o1,...,ok) divides the carmichael value of that ring, + thus there is a primitive root of order lcm(o1,...,ok). +-} +module Synthesizer.Basic.NumberTheory where + +import qualified Synthesizer.State.Signal as SigS + +import qualified Data.Set as Set +import qualified Data.Map as Map + +import qualified Algebra.Ring as Ring +import qualified Algebra.Units as Units +import qualified Algebra.PrincipalIdealDomain as PID +import qualified Algebra.IntegralDomain as Integral +import qualified Algebra.ZeroTestable as ZeroTestable + +import qualified Number.ResidueClass.Check as RC +import Number.ResidueClass.Check ((/:), ) + +import qualified Number.FixedPoint as FP +import Data.Bits (Bits, (.&.), (.|.), shiftR, ) + +import qualified Data.List.HT as ListHT +import Data.List (unfoldr, mapAccumL, genericDrop, genericSplitAt, ) +import Data.Tuple.HT (mapFst, mapSnd, mapPair, swap, ) +import Data.Maybe.HT (toMaybe, ) + +import Test.QuickCheck (Arbitrary(arbitrary), ) + +import NumericPrelude.Numeric +import NumericPrelude.Base + + +{- | +The first pair member in @powerOfTwoFactors n@ +is the maximum factor of @n@, that is a power of two. +-} +powerOfTwoFactors :: + (Bits a, Integral.C a) => a -> (a, a) +powerOfTwoFactors n = + let powerOfTwo = n .&. (-n) + in (powerOfTwo, div n powerOfTwo) + + +{- | +List all factorizations of an odd number +where the first factor is at most the second factor +and the first factors are in descending order. +-} +fermatFactors :: Integer -> [(Integer,Integer)] +fermatFactors n = + let root = FP.sqrt 1 n + in map (\(a,b) -> (b-a,b+a)) $ + mergeAndFilter + (zip (scanl (+) n [1,3..]) [0 .. div (n-1) 2]) + (zip (scanl (+) (root*root) $ iterate (2+) (2*root+1)) [root..]) + +mergeAndFilter :: (Ord a) => [(a,b)] -> [(a,c)] -> [(b,c)] +mergeAndFilter ((a0,b):a0s) ((a1,c):a1s) = + case compare a0 a1 of + LT -> mergeAndFilter a0s ((a1,c):a1s) + GT -> mergeAndFilter ((a0,b):a0s) a1s + EQ -> (b,c) : mergeAndFilter a0s a1s +mergeAndFilter _ _ = [] + + + +{- | +Argument must be a prime. +Usage of Set for efficient filtering of candidates seems to be overkill, +since the multiplicative generator seems to be small in most cases, +i.e. 2 or 3. +-} +multiplicativeGenerator :: Integer -> Integer +multiplicativeGenerator p = + let search candidates = + case Set.minView candidates of + Nothing -> error $ show p ++ " is not an odd prime" + Just (x,rest) -> + case orbitSet $ orbit p x of + new -> + -- fromIntegral (Set.size new) == p-2 + if new == Set.fromList [1..p-1] + then x + else search (Set.difference rest new) + in search (Set.fromList [2..p-1]) + + +newtype Order = Order {getOrder :: Integer} + deriving (Show, Eq, Ord) + +instance Arbitrary Order where + arbitrary = fmap (Order . (1+) . abs) arbitrary + +instance Enum Order where + succ (Order n) = Order (n+1) + pred (Order n) = Order (n-1) + fromEnum (Order n) = fromEnum n + toEnum n = Order (toEnum n) + enumFrom (Order from) = + map Order $ enumFrom from + enumFromThen (Order from) (Order thn) = + map Order $ enumFromThen from thn + enumFromTo (Order from) (Order to) = + map Order $ enumFromTo from to + enumFromThenTo (Order from) (Order thn) (Order to) = + map Order $ enumFromThenTo from thn to + +countOrder :: [a] -> Order +countOrder = foldl (\o _ -> succ o) (Order 0) + +dividesOrder :: Order -> Order -> Bool +dividesOrder (Order k) (Order n) = + divides k n + + +-- class Integral.C a => PrimitiveRoot a where +class PID.C a => PrimitiveRoot a where + primitiveRootCandidates :: a -> [a] + maximumOrderOfPrimitiveRootsOfUnity :: a -> Order + +instance PrimitiveRoot Integer where + primitiveRootCandidates modu = [1 .. modu-1] + maximumOrderOfPrimitiveRootsOfUnity = + maximumOrderOfPrimitiveRootsOfUnityInteger + +{- +For 'ordersOfPrimitiveRootsOfUnityInteger' +and the connection to Euler's totient function +we have chosen to have + +> primitiveRootsOfUnity m 1 == [1]. +-} +primitiveRootsOfUnity :: + (PrimitiveRoot a, Eq a) => a -> Order -> [a] +primitiveRootsOfUnity = + primitiveRootsOfUnityPower + +{- +Verifying that a ring has no primitive root of the wanted order +takes a long time if we do it by exhaustive search. +In the case of a=Integer we could first check, +whether the considered residue ring has a primitive root of wanted order +using the Carmichael function. +We could certainly count the number of primitive roots +and stop searching if we reach that count. +-} +primitiveRootsOfUnityPower :: + (PrimitiveRoot a, Eq a) => a -> Order -> [a] +primitiveRootsOfUnityPower modu (Order order) = + let greatDivisors = map (div order) $ uniquePrimeFactors order + in filter + (\n -> + let pow y = RC.representative $ (n /: modu) ^ y + in PID.coprime n modu + && + pow order == one + && + all (\y -> pow y /= one) greatDivisors) $ + primitiveRootCandidates modu + +primitiveRootsOfUnityNaive :: + (PrimitiveRoot a, Eq a) => a -> Order -> [a] +primitiveRootsOfUnityNaive _ (Order 0) = [] +primitiveRootsOfUnityNaive modu (Order expo) = + filter + (\n -> + let (prefix,end:_) = + genericSplitAt (expo-1) $ SigS.toList $ orbit modu n + in all (1/=) prefix && end==1) $ + primitiveRootCandidates modu + +orbitSet :: Ord a => SigS.T a -> Set.Set a +orbitSet list = + SigS.foldR + (\new cont seen -> + if Set.member new seen + then seen + else cont (Set.insert new seen)) + id list Set.empty + +orbit :: (Integral.C a) => a -> a -> SigS.T a +orbit p x = SigS.iterate (\y -> mod (x*y) p) x + + +{- | +Does not emit values in ascending order +and may return duplicates (e.g. primitiveRootsOfUnityFullOrbit 70000 10). +I hoped it would be faster than the other implementations +since it eliminates non-roots in large blocks. +However it seems that managing the root candidates in a Set +reduces performance significantly. + +The idea: +Start with a seed that is a unit. +Compute its orbit until a one is reached. +Since it is a unit, we always encounter a one. +We do not need to check for non-unit seeds, +since (gcd modu seed) will be a divisor of all seed powers. +In the orbit all numbers are powers of each other. +Now finding the roots is a matter of solving +a Diophantine equation of the exponents. +In one such step all powers in an orbit are classified as roots or non-roots +and thus we can remove them all from the set of root candidates +and continue with the remaining candidates. +Duplicates can occur if a seed +in a later iteration is found again as power of another seed. +-} +primitiveRootsOfUnityFullOrbit :: + (PrimitiveRoot a, Ord a) => a -> Order -> [a] +primitiveRootsOfUnityFullOrbit modu expo = + let search candidates = + flip fmap (Set.minView candidates) $ \(x,rest) -> + mapSnd (Set.difference rest . Set.fromList) $ + primitiveRootsOfOrbit modu expo x + in concat $ unfoldr search $ Set.fromList $ + -- needed for modules with many divisors + filter (PID.coprime modu) $ + primitiveRootCandidates modu + +primitiveRootsOfUnityFullOrbitTest :: + (PrimitiveRoot a, Ord a) => a -> Order -> [(a,[a])] +primitiveRootsOfUnityFullOrbitTest modu expo = + let search candidates = + flip fmap (Set.minView candidates) $ \(x,rest) -> + mapPair ((,) x, + Set.difference rest . Set.fromList) $ + primitiveRootsOfOrbit modu expo x + in unfoldr search $ Set.fromList $ + filter (PID.coprime modu) $ + primitiveRootCandidates modu + +primitiveRootsOfOrbit :: + (PrimitiveRoot a, Ord a) => a -> Order -> a -> ([a], [a]) +primitiveRootsOfOrbit modu (Order expo) x = + let orb = (1:) $ takeWhile (1/=) $ iterate (\y -> mod (x*y) modu) x + (Order orbitSize) = countOrder orb + in (if expo==0 + then [] + else + {- + size = length orb + Search for m and k with 0 div size expo = gcd size m + s = gcd size m + We have to consider for m + only the multiples of s. + Then divide both sides of the equation by s, yielding + 1 = gcd expo m' + -} + case divMod orbitSize expo of + (s,0) -> + map snd $ filter (PID.coprime expo . fst) $ + zip + [0 .. expo-1] + -- (ListHT.sieve s $ orb) + (map head $ iterate (genericDrop s) orb) + _ -> [], + orb) + + +hasPrimitiveRootOfUnityNaive :: + (PrimitiveRoot a, Ord a) => a -> Order -> Bool +hasPrimitiveRootOfUnityNaive modu expo = + any (dividesOrder expo . snd) $ + ordersOfPrimitiveRootsOfUnityTest modu + +{- +This should be a maximum both with respect to magnitude and to divisibility. +-} +maximumOrderOfPrimitiveRootsOfUnityNaive :: + (PrimitiveRoot a, Ord a) => a -> Order +maximumOrderOfPrimitiveRootsOfUnityNaive = + foldl max (Order 1) . map snd . ordersOfPrimitiveRootsOfUnityTest + +{- | +Computes a list of seeds and according maximum orders of roots of unity. +All divisors of those maximum orders are possible orders of roots of unity, too. +-} +ordersOfPrimitiveRootsOfUnityTest :: + (PrimitiveRoot a, Ord a) => a -> [(a, Order)] +ordersOfPrimitiveRootsOfUnityTest modu = + let search candidates = + flip fmap (Set.minView candidates) $ \(x,rest) -> + mapPair ((,) x, + Set.difference rest . Set.fromList) $ + orderOfOrbit modu x + in unfoldr search $ Set.fromList $ + filter (PID.coprime modu) $ + primitiveRootCandidates modu + +{- | +modu and x must be coprime. +If they are not, +then none of the numbers in the orbit is a root of unity +and the function enters an infinite loop. +-} +orderOfOrbit :: + (PrimitiveRoot a, Ord a) => a -> a -> (Order, [a]) +orderOfOrbit modu x = + let cyc = takeWhile (one/=) $ SigS.toList $ orbit modu x + in (succ $ countOrder cyc, cyc) + + +{- +This test speeds up 'hasPrimitiveRootOfUnityNaive' considerably +by considering the prime factors of modu. +If modu is a prime, then the ring has a multiplicative generator, +i.e. a primitive root of unity of order modu-1. +That is, all primitive roots of unity are of an order that divides modu-1. +It seems that if modu is a power of a prime, +then the according ring has also multiplicative generator. +I think this is the reason for generalising the Rader transform +to signals of prime power length. +-} +hasPrimitiveRootOfUnityInteger :: + Integer -> Order -> Bool +hasPrimitiveRootOfUnityInteger modu expo = + dividesOrder expo $ + maximumOrderOfPrimitiveRootsOfUnityInteger modu + +{- +Carmichael theorem: +If the integer residue rings with coprime moduli m0 and m1 +have primitive roots of maximum order o0 and o1, respectively, +then the integer ring with modulus m0*m1 +has maximum order (lcm o0 o1). +-} + +{- +This is the Carmichael function. +OEIS-A002322 +-} +maximumOrderOfPrimitiveRootsOfUnityInteger :: + Integer -> Order +maximumOrderOfPrimitiveRootsOfUnityInteger = + Order . + lcmMulti . + map + (\(e,p) -> + if p == 2 && e > 2 + then p^(e-2) + else p^(e-1) * (p-1)) . + map (mapFst fromIntegral) . + primeFactors + + +{- +The sum of the sub-lists should equal the Euler totient function values +(OEIS-A000010). +-} +ordersOfPrimitiveRootsOfUnityInteger :: [[Int]] +ordersOfPrimitiveRootsOfUnityInteger = + flip map [1..] $ \modu -> + let maxOrder = maximumOrderOfPrimitiveRootsOfUnity (modu::Integer) + in map (length . primitiveRootsOfUnityPower modu) $ +-- filter (flip divides maxOrder) $ + [Order 1 .. maxOrder] + +ordersOfRootsOfUnityInteger :: [[Int]] +ordersOfRootsOfUnityInteger = + flip map [1..] $ \modu -> + map (length . rootsOfUnityPower (modu::Integer)) $ + [Order 1 ..] +{- +mapM_ print $ map (\n -> (n, maximumOrderOfPrimitiveRootsOfUnityInteger (fromIntegral n), take 30 $ ordersOfRootsOfUnityInteger !! (n-1))) [2..30] + +mapM_ print $ map (\n -> (n, maximumOrderOfPrimitiveRootsOfUnityInteger (fromIntegral n), let row = ordersOfRootsOfUnityInteger !! (n-1) in map (row!!) $ map pred $ take 10 $ iterate (2*) 1)) [2..30] +-} + +ordersOfRootsOfUnityIntegerCondensed :: [[Int]] +ordersOfRootsOfUnityIntegerCondensed = + flip map [1..] $ \modu -> + let maxOrder = maximumOrderOfPrimitiveRootsOfUnity (modu::Integer) + in map (length . rootsOfUnityPower modu) $ +-- filter (flip divides maxOrder) $ + [Order 1 .. maxOrder] + +rootsOfUnityPower :: + (PrimitiveRoot a, Eq a) => a -> Order -> [a] +rootsOfUnityPower modu (Order expo) = + filter + (\n -> + PID.coprime n modu + && + RC.representative ((n /: modu) ^ expo) == one) $ + primitiveRootCandidates modu + +{- +Corollary from the Carmichael function properties: +If in Z_n there is a primitive root r of unity of order o, +then for every Z_{m \cdot n} there is also a primitive root of unity +with the same order. + +Collary: +If in Z_n1 you have a primitive root of order o1, +and so on for Z_{n_k} and order ok, +then Z_{lcm(n1,...,nk)} has primitive roots for every of the order o1,...,on. + +Conjecture: +If Z_n has a total number of m primitive roots of unity of order o, +then Z_{k*n} has at least m primitive roots of unity of order o. + +Conjecture: +Precondition: In Z_n there is a primitive root of prime order o. +Claims: +a) There are natural numbers m and k with n = m*(k*o+1) or n = m*o. +b) The smallest such n is of the form k*o+1 with k>1. + +Counterexample for a) and non-prime o: o=15, n=77 +Counterexample for b) and non-prime o: + o=20, n=25; o=27, n=81; o=54, n=81; o=55, n=121 + +Corollary from definition of Carmichael function: +For n>1, Z_{2^{n+2}} has primitive root of unity of order 2^n. +-} + +{- | +Given an order find integer residue rings +where roots of unity of this order exist. +The way they are constructed also warrants, +that 'order' is a unit (i.e. invertible) in those rings. + +The list is not exhaustive +but computes suggestions quickly. +The first found modulus seems to be smallest one that exist. +However, the first modulus is not the smallest one +among the ones that only have the wanted primitive root, +but where 'order' is not necessarily a unit. +E.g. + +> ringsWithPrimitiveRootOfUnityAndUnit 840 == 2521 : 3361 : ... + +but the smallest modulus is 1763. + +Most of the numbers are primes. +For these the recursion property of the Carmichael function +immediately yields that there are primitive roots of unity of order 'order'. +-} +ringsWithPrimitiveRootOfUnityAndUnit :: Order -> [Integer] +ringsWithPrimitiveRootOfUnityAndUnit order@(Order k) = + filter (flip hasPrimitiveRootOfUnityInteger order) $ + iterate (k+) 1 + + +ringsWithPrimitiveRootsOfUnityAndUnitsNaive :: [Order] -> [Integer] -> [Integer] +ringsWithPrimitiveRootsOfUnityAndUnitsNaive rootOrders units = + filter + (\n -> + all (hasPrimitiveRootOfUnityInteger n) rootOrders && + all (PID.coprime n) units) + [1..] + + +{- +It would be nice to have the Omega monad here +in order to enumerate all rings. +-} +ringWithPrimitiveRootsOfUnityAndUnits :: [Order] -> [Integer] -> Integer +ringWithPrimitiveRootsOfUnityAndUnits rootOrders units = + let p = lcmMulti units + in lcmMulti $ + map (head . filter (PID.coprime p) . + ringsWithPrimitiveRootOfUnityAndUnit) $ + rootOrders + +{- +Search for an appriopriate modulus +using the recursive definition of the Carmichael function. +We chose the prime factors of the Carmichael function argument +such that we get at least the prime factors in the function value that we need. + +The modulus constructed this way is usually not the smallest possible +and it also does not provide that 'n' is a unit in the residue ring. +The simpler function 'ringsWithPrimitiveRootOfUnityAndUnit' +will usually produce a smaller modulus. +-} +ringWithPrimitiveRootsOfUnity :: Order -> Integer +ringWithPrimitiveRootsOfUnity (Order n) = + case n of + 0 -> 2 + _ -> + product $ map (uncurry ringPower) $ snd $ + mapAccumL + (\factors (e,p) -> + if Map.findWithDefault 0 p factors >= e + then (factors, (0,p)) + else + if p==2 + then + (factors, + case e of + 0 -> (0,2) + 1 -> (1,3) + 2 -> (1,5) + _ -> (e+2, 2)) + else + (Map.unionWith max factors $ + Map.fromList $ map swap $ primeFactors $ p-1, + (e+1, p))) + Map.empty $ + reverse $ primeFactors $ lcmMulti $ + n : map (subtract 1) (partialPrimes n) + +lcmMulti :: (PID.C a) => [a] -> a +lcmMulti = foldl lcm one + + +{- | +List all numbers that only contain prime factors 2 and 3 in ascending order. +OEIS:A003586 +-} +numbers3Smooth :: [Integer] +numbers3Smooth = + foldr + (\(x0:x1:xs) ys -> x0 : x1 : ListHT.mergeBy (<=) xs ys) + (error "numbers3Smooth: infinite list should not have an end") $ + iterate (map (3*)) $ + iterate (2*) 1 + +numbers3SmoothAlt :: [Integer] +numbers3SmoothAlt = + unfoldr + (fmap (\(m,rest) -> (m, Set.union rest $ Set.fromAscList [2*m,3*m])) . + Set.minView) $ + Set.singleton 1 + +{- +OEIS:A051037 +-} +numbers5Smooth :: [Integer] +numbers5Smooth = + foldr + (\(x0:x1:x2:xs) ys -> x0 : x1 : x2 : ListHT.mergeBy (<=) xs ys) + (error "numbers5Smooth: infinite list should not have an end") $ + iterate (map (5*)) $ + numbers3Smooth + +numbers5SmoothAlt :: [Integer] +numbers5SmoothAlt = + unfoldr + (fmap (\(m,rest) -> (m, Set.union rest $ Set.fromAscList [2*m,3*m,5*m])) . + Set.minView) $ + Set.singleton 1 + +ceilingPowerOfTwo :: (Ring.C a, Bits a) => a -> a +ceilingPowerOfTwo 0 = 1 +ceilingPowerOfTwo n = + (1+) $ fst $ head $ + dropWhile (uncurry (/=)) $ + ListHT.mapAdjacent (,) $ + scanl (\m d -> shiftR m d .|. m) (n-1) $ + iterate (2*) 1 + +divideByMaximumPower :: + (Integral.C a, ZeroTestable.C a) => a -> a -> a +divideByMaximumPower b n = + last $ + n : unfoldr (\m -> case divMod m b of (q,r) -> toMaybe (isZero r) (q,q)) n + +divideByMaximumPowerRecursive :: + (Integral.C a, Eq a, ZeroTestable.C a) => a -> a -> a +divideByMaximumPowerRecursive b = + let recourse n = + case divMod b n of + (q,0) -> recourse q + _ -> n + in recourse + +getMaximumExponent :: + (Integral.C a, ZeroTestable.C a) => + a -> a -> (Int,a) +getMaximumExponent b n = + last $ (0,n) : + unfoldr + (\(e,m) -> + let (q,r) = divMod m b + eq = (e+1,q) + in toMaybe (isZero r) (eq,eq)) + (0,n) + +is3Smooth :: Integer -> Bool +is3Smooth = + (1==) . + divideByMaximumPower 3 . + divideByMaximumPower 2 + +is5Smooth :: Integer -> Bool +is5Smooth = + (1==) . + divideByMaximumPower 5 . + divideByMaximumPower 3 . + divideByMaximumPower 2 + +{- | +Compute the smallest composite of 2 and 3 that is as least as large as the input. +This can be interpreted as solving an integer linear programming problem with +min (\(a,b) -> a * log 2 + b * log 3) +over the domain {(a,b) : a>=0, b>=0, a * log 2 + b * log 3 >= log n} +-} +{- +Problem: We cannot just start with the ceilingPowerOfTwo +and then multiply with 3/4 until we fall below n, +since the 3/4 decreases too fast. +27/32 is closer to one, +and higher powers of 3 and 2 in the ratio make the ratio even closer to one. +-} +ceiling3Smooth :: Integer -> Integer +ceiling3Smooth n = + head $ dropWhile ( Integer +ceiling5Smooth n = + head $ dropWhile ( Integer +ceiling3SmoothNaive = + head . + dropWhile (not . is3Smooth) . + iterate (1+) + +ceiling5SmoothNaive :: Integer -> Integer +ceiling5SmoothNaive = + head . + dropWhile (not . is5Smooth) . + iterate (1+) + + +{- | +Compute all primes that occur in the course of dividing +a Fourier transform into sub-transforms. +-} +partialPrimes :: Integer -> [Integer] +partialPrimes = + let primeFactorSet = + Set.fromAscList . uniquePrimeFactors + in unfoldr + (fmap + (\(p,set) -> + (p, Set.union (primeFactorSet (p-1)) set)) . + Set.maxView) + . + primeFactorSet + +-- cf. htam:NumberTheory +uniquePrimeFactors :: + (Integral.C a, Bits a, ZeroTestable.C a, Ord a) => + a -> [a] +-- map snd . primeFactors +uniquePrimeFactors n = + let oddFactors = + foldr + (\p go m -> + let (q,r) = divMod m p + in if r==0 + then p : go (divideByMaximumPower p q) + else + if q >= p + then go m + else if m==1 then [] else m : []) + (error "uniquePrimeFactors: end of infinite list") + (iterate (2+) 3) + in case powerOfTwoFactors n of + (1,m) -> oddFactors m + (_,m) -> 2 : oddFactors m + +{- | +Prime factors and their exponents in ascending order. +-} +primeFactors :: + (PrimitiveRoot a, Ord a) => a -> [(Int, a)] +primeFactors n = + let oddFactors = + foldr + (\p go m -> + let (q0,r) = divMod m p + in if r==0 + then + case getMaximumExponent p q0 of + (e,q1) -> (e+1,p) : go q1 + else + if q0 >= p + then go m + else if m==1 then [] else (1,m) : []) + (const []) + (filter (not . Units.isUnit) $ + primitiveRootCandidates n) + in case getMaximumExponent 2 n of + (0,m) -> oddFactors m + (e,m) -> (e,2) : oddFactors m + +{- +cf. htam:NumberTheory + +Shall this be moved to NumericPrelude? + +It should be replaced by a more sophisticated prime test. +-} +isPrime :: Integer -> Bool +isPrime n = + case primeFactors n of + [] -> False + (e,m):_ -> e==1 && m==n + +{- | +Find lengths of signals that require many interim Rader transforms +and end with the given length. + +raderWorstCases 2 = OEIS-A061092 +raderWorstCases 5 = tail OEIS-A059411 + +Smallest raderWorstCase numbers are 2,5,13,17,19,31,37,41,43,61,... +This matches the definition of OEIS-A061303. +-} +raderWorstCases :: Integer -> [Integer] +raderWorstCases = + iterate + (\n -> + head $ dropWhile (not . isPrime) $ + tail $ iterate (n+) 1) + +{- | +This is usually faster than 'fastFourierRing' +since it does not need to factor large numbers. +However, the generated modulus is usually much greater. +-} +{- +I see the following opportunities for optimization: + +1. Speedup 'fastFourierRing' by + faster primality test (Miller-Rabin) and + faster prime factorization (Pollard-Rho-method). + These are also needed for + maximumOrderOfPrimitiveRootsOfUnityInteger + that is used by Fourier.Element.primitiveRoot + in order to compute a root with maximum order. + +2. Reduce the moduli produced by 'fastFourierRingAlt' + by merging some orders that are passed to + ringWithPrimitiveRootsOfUnityAndUnits, + such that an LCM of a group of orders can still be handled. + This is a kind of knapsack problem. + Maybe we could collect the factors in a way + such that (lcm orderGroup + 1) is prime. + +3. Avoid to compute factorizations of numbers + where we already know the factors + or where we do not need the factors at all. + Use the factors returned by partialPrimes + in order to compute a prime factorization + of lcmMulti (map pred (partialPrimes n)). + Call this (product ps). + Now search for rings of moduli (1 + k * product ps), + where there are (small) primitive roots of order (product ps). + We only need to check whether there are small numbers + such as 2, 3, 5, 6, 7 that have a (product ps)-th power that is 1, + using fast exponentiation. + If there is a power being 1 then check for primitivity + by computing (k * product ps / p)-th powers + for all prime factors p of (k * product ps). + If there is no primitive root <= 7, + there may still be a primitive root of wanted order, + but it is then cheaper to seek for larger moduli. + + If we finally have a nice modulus + it is still stupid to factorize (modulus-1) + and search for a primitive root + in each invocation of Fourier.Element.primitiveRoot. + We could define a special datatype analogously to ResidueClass, + that stores the primitive root and its order + additional to the ResidueClass modulus. +-} +fastFourierRingAlt :: Int -> Integer +fastFourierRingAlt n = + case n of + 0 -> 2 + 1 -> 2 + _ -> + let ni = fromIntegral n + ps = filter (>1) (map (subtract 1) (partialPrimes ni)) + in ringWithPrimitiveRootsOfUnityAndUnits (map Order $ ni : ps) ps + +{- | +Determine an integer residue ring +in which a Fast Fourier transform of size n can be performed. +It must contain certain primitive roots. +If we choose a non-primitive root, +then some off-diagonal values in F^-1·F are non-zero. +-} +{- +When we need roots of orders o1,...,ok and according inverse elements +we can also ask for a ring, where there is a root of order lcm(o1,...,ok). +The answer to both questions is the same set of rings. +This can be proven using the statement, +that the order of any primitive root +divides the carmichael value of the modulus. + +Since ringWithPrimitiveRootsOfUnityAndUnits +multiplies the moduli of rings for o1,...,ok, +it will produce large moduli. +-} +fastFourierRing :: Int -> Integer +fastFourierRing n = + case n of + 0 -> 2 + 1 -> 2 + _ -> + let ni = fromIntegral n + in {- + We cannot use ringsWithPrimitiveRootOfUnityAndUnit + since for 359 we already get an Int overflow. + For 719, 1439, 2879 we also get a very large value. + -} + head $ filter isPrime $ + (\order -> iterate (order +) 1) $ + lcmMulti $ + ni : map (subtract 1) (partialPrimes ni) diff --git a/private/Synthesizer/Generic/Permutation.hs b/private/Synthesizer/Generic/Permutation.hs new file mode 100644 --- /dev/null +++ b/private/Synthesizer/Generic/Permutation.hs @@ -0,0 +1,151 @@ +{- | +Permutations of signals as needed for Fast Fourier transforms. +Most functions are independent of the Signal framework. +We could move them as well to Synthesizer.Basic. +-} +module Synthesizer.Generic.Permutation where + +import qualified Synthesizer.Basic.NumberTheory as NumberTheory + +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.State.Signal as SigS + +import qualified Data.StorableVector.ST.Strict as SVST +import qualified Data.StorableVector as SV + +import qualified Algebra.PrincipalIdealDomain as PID + + + +type T = SV.Vector Int + +apply :: + (SigG.Transform sig y) => + T -> sig y -> sig y +apply p xs = + SigG.takeStateMatch xs $ + SigS.map (SigG.index xs) $ + SigS.fromStrictStorableSignal p + + +size :: T -> Int +size = SV.length + + +{- | +> inverse (transposition n m) = transposition m n +-} +transposition :: + Int -> Int -> T +transposition n m = + fst $ SV.unfoldrN (n*m) + (\(i,j,k0) -> Just (i, + case pred k0 of + 0 -> let j1 = j+1 in (j1, j1, m) + k1 -> (i+n, j, k1))) + (0,0,m) + + +{- +In general the inverse of a skewGrid +does not look like even a generalized skewGrid. +E.g. @inverse $ skewGrid 3 4@. +-} +skewGrid :: + Int -> Int -> T +skewGrid n m = + let len = n*m + in fst $ SV.unfoldrN len + (\(i0,k0) -> Just (i0, + let k1 = pred k0 + i1 = i0+n + in if k1==0 + then (mod (i1+m) len, m) + else (mod i1 len, k1))) + (0,m) + +{- | +> inverse (skewGrid n m) == skewGridInv n m + +In general the inverse of a skewGrid +cannot be expressed like skewGrid or skewGridCRT. +E.g. @inverse $ skewGrid 3 4@. +-} +skewGridInv :: + Int -> Int -> T +skewGridInv n m = + SV.pack $ + map (\k -> + let Just (i,j) = PID.diophantine k n m + in mod i m + mod j n * m) $ + take (n*m) $ iterate (1+) 0 + +skewGridCRT :: + Int -> Int -> T +skewGridCRT n m = + let len = n*m + (ni,mi) = snd $ PID.extendedGCD n m + in fst $ SV.unfoldrN len + (\(i0,k0) -> Just (i0, + let k1 = pred k0 + i1 = i0+ni*n + in if k1==0 + then (mod (i1+mi*m) len, m) + else (mod i1 len, k1))) + (0,m) + +skewGridCRTInv :: + Int -> Int -> T +skewGridCRTInv n m = + fst $ SV.packN (n*m) $ + map (\k -> mod k m + mod k n * m) $ + iterate (1+) 0 + + +{- | +Beware of 0-based indices stored in the result vector. +-} +multiplicative :: Int -> T +multiplicative ni = + let n = fromIntegral ni + gen = NumberTheory.multiplicativeGenerator n + in {- + Since 'gen' is usually 2 or 3, + the error should occur really only for huge signals. + -} + if gen * n > fromIntegral (maxBound :: Int) + then error "signal too long for Int indexing" + else fst $ SV.unfoldrN (ni-1) + (\x -> Just (x-1, mod (fromInteger gen * x) ni)) 1 + +{- | +We only need to compute the inverse permutation explicitly, +because not all signal structures support write to arbitrary indices, +thus Generic.Write does not support it. +For strict StorableVector it would be more efficient +to build the vector directly. + +It holds: + +> inverse . inverse == id +-} +inverse :: T -> T +inverse perm = + SVST.runSTVector + (do inv <- SVST.new_ (SV.length perm) + SigS.sequence_ $ + SigS.zipWith (SVST.write inv) + (SigS.fromStrictStorableSignal perm) + (SigS.iterate (1+) 0) + return inv) + +reverse :: T -> T +reverse perm = + fst $ SV.unfoldrN (SV.length perm) + (\mn -> Just $ + case mn of + Nothing -> (SV.head perm, Just $ SV.length perm) + Just n -> + let n1 = n-1 + in (SV.index perm n1, Just n1)) + Nothing diff --git a/speedtest/SpeedTest.hs b/speedtest/SpeedTest.hs --- a/speedtest/SpeedTest.hs +++ b/speedtest/SpeedTest.hs @@ -11,8 +11,7 @@ import qualified Data.Binary.Put as Bin import Foreign (Int16, Ptr, alloca, allocaBytes, poke, pokeElemOff, sizeOf) -import System.IO (openBinaryFile, IOMode(WriteMode), hClose, Handle, hPutBuf) -import Control.Exception (bracket) +import System.IO (withBinaryFile, IOMode(WriteMode), Handle, hPutBuf) import qualified Algebra.Transcendental as Trans import qualified Algebra.RealField as RealField @@ -172,7 +171,7 @@ writeSignalMonoPoke :: FilePath -> [Int16] -> IO () writeSignalMonoPoke fileName signal = - bracket (openBinaryFile fileName WriteMode) hClose $ + withBinaryFile fileName WriteMode $ \h -> alloca $ \p -> mapM_ (putInt h p) signal @@ -190,7 +189,7 @@ writeSignalMonoBlock :: FilePath -> [Int16] -> IO () writeSignalMonoBlock fileName signal = - bracket (openBinaryFile fileName WriteMode) hClose $ + withBinaryFile fileName WriteMode $ \h -> let blocks = sliceVertical maxBlockSize signal in allocaBytes (int16size * maxBlockSize) $ \p -> mapM_ (putIntBlock h p) blocks @@ -214,7 +213,7 @@ writeZeroBlocks :: FilePath -> Int -> IO () writeZeroBlocks fileName len = - bracket (openBinaryFile fileName WriteMode) hClose $ + withBinaryFile fileName WriteMode $ \h -> allocaBytes (int16size * maxBlockSize) $ \p -> do mapM_ (\off -> pokeElemOff p off (P98.fromInteger 0 :: Int16)) diff --git a/speedtest/SpeedTestExp.hs b/speedtest/SpeedTestExp.hs --- a/speedtest/SpeedTestExp.hs +++ b/speedtest/SpeedTestExp.hs @@ -9,14 +9,13 @@ import qualified Data.ByteString.Lazy as B import qualified Data.Binary.Put as Bin -import Data.Array.IO (IOUArray, newArray_, castIOUArray, hPutArray, writeArray) +import Data.Array.IO (IOUArray, newArray_, hPutArray, writeArray) +import Data.Array.Unsafe (castIOUArray) import Data.Word(Word8) --- we could also use withBinaryFile -import System.IO (openBinaryFile, hClose, hPutBuf, IOMode(WriteMode)) +import System.IO (withBinaryFile, hPutBuf, IOMode(WriteMode)) import Foreign (Int16, pokeElemOff, allocaBytes) -import Control.Exception (bracket) import Control.Monad (zipWithM_) import GHC.Float (double2Int) @@ -60,7 +59,7 @@ writeSignal :: FilePath -> Int -> [Double] -> IO () writeSignal name num signal = - bracket (openBinaryFile name WriteMode) hClose $ \h -> + withBinaryFile name WriteMode $ \h -> allocaBytes (2*num) $ \buf -> zipWithM_ (pokeElemOff buf) [0..(num-1)] @@ -69,7 +68,7 @@ writeExponentialList :: FilePath -> Int -> Double -> Double -> IO () writeExponentialList name num hl y0 = - bracket (openBinaryFile name WriteMode) hClose $ \h -> + withBinaryFile name WriteMode $ \h -> allocaBytes (2*num) $ \buf -> zipWithM_ (pokeElemOff buf) [0..(num-1)] @@ -79,7 +78,7 @@ writeExponential :: FilePath -> Int -> Double -> Double -> IO () writeExponential name num hl y0 = - bracket (openBinaryFile name WriteMode) hClose $ \h -> + withBinaryFile name WriteMode $ \h -> allocaBytes (2*num) $ \buf -> {- let k = 0.5**(1/hl) @@ -103,7 +102,7 @@ writeExponentialIOUArray :: FilePath -> Int -> Double -> Double -> IO () writeExponentialIOUArray name num hl y0 = - bracket (openBinaryFile name WriteMode) hClose $ \h -> + withBinaryFile name WriteMode $ \h -> newArray_ (0,2*num-1) >>= \arr -> let k = 0.5**(1/hl) loop i y = @@ -118,7 +117,7 @@ writeExponentialStorableVector :: FilePath -> Int -> Double -> Double -> IO () writeExponentialStorableVector name num hl y0 = - bracket (openBinaryFile name WriteMode) hClose $ \h -> + withBinaryFile name WriteMode $ \h -> let k = 0.5**(1/hl) (fp, _offset, _size) = VB.toForeignPtr $ fst $ diff --git a/src/Synthesizer/Basic/NumberTheory.hs b/src/Synthesizer/Basic/NumberTheory.hs deleted file mode 100644 --- a/src/Synthesizer/Basic/NumberTheory.hs +++ /dev/null @@ -1,896 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -{- -Some of these functions might be moved to NumericPrelude. - -Wikipedia: (primitive) roots of unity modulo n - (primitive) roots must be units and all units are (primitive) roots - maximum possible order for primitive roots - Carmichael - all possible orders: divisor of Carmichael (proof? statement already in Carmichael-function-article) - sum of primitive roots that vanishes - order of primitive root is a divisor of each possible exponent - proof with GCD and diophantine in exponent - check for primitive root: fast exponentiation, - primitivity: check exponents that are prime divisors - how to find a primitive root: just try - sum of powers of a primitive root is zero - number of primitive roots of given order - g(n,k) > 0 if k|lambda(n) - g(n,k) = 0 else - g(n,1) = 1 - g(4,2) = 1 - g(2^n,2) = 3 for n>=3 ((-1) is always a square root of 1) - g(2^n,2^k) = 2^k for k>=2 && k=3 and n in OEIS:A033948 - sum(g(n,k), k\in\N) = phi(n) - There are only a few patterns that occur as rows of g, - but a row of g (i.e. g(n)) does functionally depend on - either lambda(n) nor phi(n) - lambda(14) = 6 nozeros(g(14)) = [1,1,2,2] (f ~ [1,2,3,6]) - lambda(21) = 6 nozeros(g(21)) = [1,3,2,6] (f ~ [1,4,3,12]) - phi(13) = 12 nozeros(g(13)) = [1,1,2,2,2,4] (f ~ [1,2,3,4,6,12]) - phi(21) = 12 nozeros(g(21)) = [1,3,2,6] (f ~ [1,4,3,12]) - However length(nozeros(f(n))) = numberofdivisors(lambda(n)) - numberofdivisors=A000005 - number of roots of given order - easier to compute - k|m => f(n,k) | f(n,m) - g(n,k) = f(n,k) - sum(f(n,d), d|k and k/d prime) + ... - inclusion-exclusion-principle - better to write the other round: - f(n,k) = sum(g(n,d), d|k) - this is Dirichlet convolution - RUNM says f(n,k) is multiplicative - list it in multiplicative function - f(n,1) = 1 for n>=2 - f(n,lambda(n)) = phi(n) - f(n,a·b) = f(n,a)·f(n,b) if a and b are coprime (conjecture) - f(n,lcm(a,b)) = lcm(f(n,a),f(n,b)) (conjecture) - If this conjecture is true, then we only need to know f(n,p^i). - The following conjecture is wrong: - for prime p it is f(n,p^i) = gcd(lambda(n),p^i) - counterexamples - f(8,2) = 4, lambda(8)=2 - f(63,3) = 9, lambda(63)=6 - f(275,5) = 25, lambda(275)=20 - f(1247,7) = 49, lambda(1247)=84 - It seems to be: - for prime p it is f(n,p^i) = p^j for some j - How to find a modulus where there is a primitive root of order o? - just try numbers from the sequence o+1, 2*o+1, 3*o+1 - Because of [[Dirichlet's theorem on arithmetic progressions]] - you will somewhen find a prime p, - and its Carmichael value is p-1, which is a multiple of o. - In this ring even 'o' is a unit. - How to find a modulus where there are primitive roots of orders o1,..,ok? - Just search for a modulus with roots of order lcm(o1,...,ok). - The smallest such modulus should also be the smallest one - that has primitive roots of orders o1,..,ok? - Proof: If a ring has primitive roots of orders o1,..,ok - then all orders divide the carmichael value of that ring, - thus lcm(o1,...,ok) divides the carmichael value of that ring, - thus there is a primitive root of order lcm(o1,...,ok). --} -module Synthesizer.Basic.NumberTheory where - -import qualified Synthesizer.State.Signal as SigS - -import qualified Data.Set as Set -import qualified Data.Map as Map - -import qualified Algebra.Ring as Ring -import qualified Algebra.Units as Units -import qualified Algebra.PrincipalIdealDomain as PID -import qualified Algebra.IntegralDomain as Integral -import qualified Algebra.ZeroTestable as ZeroTestable - -import qualified Number.ResidueClass.Check as RC -import Number.ResidueClass.Check ((/:), ) - -import qualified Number.FixedPoint as FP -import Data.Bits (Bits, (.&.), (.|.), shiftR, ) - -import qualified Data.List.HT as ListHT -import Data.List (unfoldr, mapAccumL, genericDrop, genericSplitAt, ) -import Data.Tuple.HT (mapFst, mapSnd, mapPair, swap, ) -import Data.Maybe.HT (toMaybe, ) - -import Test.QuickCheck (Arbitrary(arbitrary), ) - -import NumericPrelude.Numeric -import NumericPrelude.Base - - -{- | -The first pair member in @powerOfTwoFactors n@ -is the maximum factor of @n@, that is a power of two. --} -powerOfTwoFactors :: - (Bits a, Integral.C a) => a -> (a, a) -powerOfTwoFactors n = - let powerOfTwo = n .&. (-n) - in (powerOfTwo, div n powerOfTwo) - - -{- | -List all factorizations of an odd number -where the first factor is at most the second factor -and the first factors are in descending order. --} -fermatFactors :: Integer -> [(Integer,Integer)] -fermatFactors n = - let root = FP.sqrt 1 n - in map (\(a,b) -> (b-a,b+a)) $ - mergeAndFilter - (zip (scanl (+) n [1,3..]) [0 .. div (n-1) 2]) - (zip (scanl (+) (root*root) $ iterate (2+) (2*root+1)) [root..]) - -mergeAndFilter :: (Ord a) => [(a,b)] -> [(a,c)] -> [(b,c)] -mergeAndFilter ((a0,b):a0s) ((a1,c):a1s) = - case compare a0 a1 of - LT -> mergeAndFilter a0s ((a1,c):a1s) - GT -> mergeAndFilter ((a0,b):a0s) a1s - EQ -> (b,c) : mergeAndFilter a0s a1s -mergeAndFilter _ _ = [] - - - -{- | -Argument must be a prime. -Usage of Set for efficient filtering of candidates seems to be overkill, -since the multiplicative generator seems to be small in most cases, -i.e. 2 or 3. --} -multiplicativeGenerator :: Integer -> Integer -multiplicativeGenerator p = - let search candidates = - case Set.minView candidates of - Nothing -> error $ show p ++ " is not an odd prime" - Just (x,rest) -> - case orbitSet $ orbit p x of - new -> - -- fromIntegral (Set.size new) == p-2 - if new == Set.fromList [1..p-1] - then x - else search (Set.difference rest new) - in search (Set.fromList [2..p-1]) - - -newtype Order = Order {getOrder :: Integer} - deriving (Show, Eq, Ord) - -instance Arbitrary Order where - arbitrary = fmap (Order . (1+) . abs) arbitrary - -instance Enum Order where - succ (Order n) = Order (n+1) - pred (Order n) = Order (n-1) - fromEnum (Order n) = fromEnum n - toEnum n = Order (toEnum n) - enumFrom (Order from) = - map Order $ enumFrom from - enumFromThen (Order from) (Order thn) = - map Order $ enumFromThen from thn - enumFromTo (Order from) (Order to) = - map Order $ enumFromTo from to - enumFromThenTo (Order from) (Order thn) (Order to) = - map Order $ enumFromThenTo from thn to - -countOrder :: [a] -> Order -countOrder = foldl (\o _ -> succ o) (Order 0) - -dividesOrder :: Order -> Order -> Bool -dividesOrder (Order k) (Order n) = - divides k n - - --- class Integral.C a => PrimitiveRoot a where -class PID.C a => PrimitiveRoot a where - primitiveRootCandidates :: a -> [a] - maximumOrderOfPrimitiveRootsOfUnity :: a -> Order - -instance PrimitiveRoot Integer where - primitiveRootCandidates modu = [1 .. modu-1] - maximumOrderOfPrimitiveRootsOfUnity = - maximumOrderOfPrimitiveRootsOfUnityInteger - -{- -For 'ordersOfPrimitiveRootsOfUnityInteger' -and the connection to Euler's totient function -we have chosen to have - -> primitiveRootsOfUnity m 1 == [1]. --} -primitiveRootsOfUnity :: - (PrimitiveRoot a, Eq a) => a -> Order -> [a] -primitiveRootsOfUnity = - primitiveRootsOfUnityPower - -{- -Verifying that a ring has no primitive root of the wanted order -takes a long time if we do it by exhaustive search. -In the case of a=Integer we could first check, -whether the considered residue ring has a primitive root of wanted order -using the Carmichael function. -We could certainly count the number of primitive roots -and stop searching if we reach that count. --} -primitiveRootsOfUnityPower :: - (PrimitiveRoot a, Eq a) => a -> Order -> [a] -primitiveRootsOfUnityPower modu (Order order) = - let greatDivisors = map (div order) $ uniquePrimeFactors order - in filter - (\n -> - let pow y = RC.representative $ (n /: modu) ^ y - in PID.coprime n modu - && - pow order == one - && - all (\y -> pow y /= one) greatDivisors) $ - primitiveRootCandidates modu - -primitiveRootsOfUnityNaive :: - (PrimitiveRoot a, Eq a) => a -> Order -> [a] -primitiveRootsOfUnityNaive _ (Order 0) = [] -primitiveRootsOfUnityNaive modu (Order expo) = - filter - (\n -> - let (prefix,end:_) = - genericSplitAt (expo-1) $ SigS.toList $ orbit modu n - in all (1/=) prefix && end==1) $ - primitiveRootCandidates modu - -orbitSet :: Ord a => SigS.T a -> Set.Set a -orbitSet list = - SigS.foldR - (\new cont seen -> - if Set.member new seen - then seen - else cont (Set.insert new seen)) - id list Set.empty - -orbit :: (Integral.C a) => a -> a -> SigS.T a -orbit p x = SigS.iterate (\y -> mod (x*y) p) x - - -{- | -Does not emit values in ascending order -and may return duplicates (e.g. primitiveRootsOfUnityFullOrbit 70000 10). -I hoped it would be faster than the other implementations -since it eliminates non-roots in large blocks. -However it seems that managing the root candidates in a Set -reduces performance significantly. - -The idea: -Start with a seed that is a unit. -Compute its orbit until a one is reached. -Since it is a unit, we always encounter a one. -We do not need to check for non-unit seeds, -since (gcd modu seed) will be a divisor of all seed powers. -In the orbit all numbers are powers of each other. -Now finding the roots is a matter of solving -a Diophantine equation of the exponents. -In one such step all powers in an orbit are classified as roots or non-roots -and thus we can remove them all from the set of root candidates -and continue with the remaining candidates. -Duplicates can occur if a seed -in a later iteration is found again as power of another seed. --} -primitiveRootsOfUnityFullOrbit :: - (PrimitiveRoot a, Ord a) => a -> Order -> [a] -primitiveRootsOfUnityFullOrbit modu expo = - let search candidates = - flip fmap (Set.minView candidates) $ \(x,rest) -> - mapSnd (Set.difference rest . Set.fromList) $ - primitiveRootsOfOrbit modu expo x - in concat $ unfoldr search $ Set.fromList $ - -- needed for modules with many divisors - filter (PID.coprime modu) $ - primitiveRootCandidates modu - -primitiveRootsOfUnityFullOrbitTest :: - (PrimitiveRoot a, Ord a) => a -> Order -> [(a,[a])] -primitiveRootsOfUnityFullOrbitTest modu expo = - let search candidates = - flip fmap (Set.minView candidates) $ \(x,rest) -> - mapPair ((,) x, - Set.difference rest . Set.fromList) $ - primitiveRootsOfOrbit modu expo x - in unfoldr search $ Set.fromList $ - filter (PID.coprime modu) $ - primitiveRootCandidates modu - -primitiveRootsOfOrbit :: - (PrimitiveRoot a, Ord a) => a -> Order -> a -> ([a], [a]) -primitiveRootsOfOrbit modu (Order expo) x = - let orb = (1:) $ takeWhile (1/=) $ iterate (\y -> mod (x*y) modu) x - (Order orbitSize) = countOrder orb - in (if expo==0 - then [] - else - {- - size = length orb - Search for m and k with 0 div size expo = gcd size m - s = gcd size m - We have to consider for m - only the multiples of s. - Then divide both sides of the equation by s, yielding - 1 = gcd expo m' - -} - case divMod orbitSize expo of - (s,0) -> - map snd $ filter (PID.coprime expo . fst) $ - zip - [0 .. expo-1] - -- (ListHT.sieve s $ orb) - (map head $ iterate (genericDrop s) orb) - _ -> [], - orb) - - -hasPrimitiveRootOfUnityNaive :: - (PrimitiveRoot a, Ord a) => a -> Order -> Bool -hasPrimitiveRootOfUnityNaive modu expo = - any (dividesOrder expo . snd) $ - ordersOfPrimitiveRootsOfUnityTest modu - -{- -This should be a maximum both with respect to magnitude and to divisibility. --} -maximumOrderOfPrimitiveRootsOfUnityNaive :: - (PrimitiveRoot a, Ord a) => a -> Order -maximumOrderOfPrimitiveRootsOfUnityNaive = - foldl max (Order 1) . map snd . ordersOfPrimitiveRootsOfUnityTest - -{- | -Computes a list of seeds and according maximum orders of roots of unity. -All divisors of those maximum orders are possible orders of roots of unity, too. --} -ordersOfPrimitiveRootsOfUnityTest :: - (PrimitiveRoot a, Ord a) => a -> [(a, Order)] -ordersOfPrimitiveRootsOfUnityTest modu = - let search candidates = - flip fmap (Set.minView candidates) $ \(x,rest) -> - mapPair ((,) x, - Set.difference rest . Set.fromList) $ - orderOfOrbit modu x - in unfoldr search $ Set.fromList $ - filter (PID.coprime modu) $ - primitiveRootCandidates modu - -{- | -modu and x must be coprime. -If they are not, -then none of the numbers in the orbit is a root of unity -and the function enters an infinite loop. --} -orderOfOrbit :: - (PrimitiveRoot a, Ord a) => a -> a -> (Order, [a]) -orderOfOrbit modu x = - let cyc = takeWhile (one/=) $ SigS.toList $ orbit modu x - in (succ $ countOrder cyc, cyc) - - -{- -This test speeds up 'hasPrimitiveRootOfUnityNaive' considerably -by considering the prime factors of modu. -If modu is a prime, then the ring has a multiplicative generator, -i.e. a primitive root of unity of order modu-1. -That is, all primitive roots of unity are of an order that divides modu-1. -It seems that if modu is a power of a prime, -then the according ring has also multiplicative generator. -I think this is the reason for generalising the Rader transform -to signals of prime power length. --} -hasPrimitiveRootOfUnityInteger :: - Integer -> Order -> Bool -hasPrimitiveRootOfUnityInteger modu expo = - dividesOrder expo $ - maximumOrderOfPrimitiveRootsOfUnityInteger modu - -{- -Carmichael theorem: -If the integer residue rings with coprime moduli m0 and m1 -have primitive roots of maximum order o0 and o1, respectively, -then the integer ring with modulus m0*m1 -has maximum order (lcm o0 o1). --} - -{- -This is the Carmichael function. -OEIS-A002322 --} -maximumOrderOfPrimitiveRootsOfUnityInteger :: - Integer -> Order -maximumOrderOfPrimitiveRootsOfUnityInteger = - Order . - lcmMulti . - map - (\(e,p) -> - if p == 2 && e > 2 - then p^(e-2) - else p^(e-1) * (p-1)) . - map (mapFst fromIntegral) . - primeFactors - - -{- -The sum of the sub-lists should equal the Euler totient function values -(OEIS-A000010). --} -ordersOfPrimitiveRootsOfUnityInteger :: [[Int]] -ordersOfPrimitiveRootsOfUnityInteger = - flip map [1..] $ \modu -> - let maxOrder = maximumOrderOfPrimitiveRootsOfUnity (modu::Integer) - in map (length . primitiveRootsOfUnityPower modu) $ --- filter (flip divides maxOrder) $ - [Order 1 .. maxOrder] - -ordersOfRootsOfUnityInteger :: [[Int]] -ordersOfRootsOfUnityInteger = - flip map [1..] $ \modu -> - map (length . rootsOfUnityPower (modu::Integer)) $ - [Order 1 ..] -{- -mapM_ print $ map (\n -> (n, maximumOrderOfPrimitiveRootsOfUnityInteger (fromIntegral n), take 30 $ ordersOfRootsOfUnityInteger !! (n-1))) [2..30] - -mapM_ print $ map (\n -> (n, maximumOrderOfPrimitiveRootsOfUnityInteger (fromIntegral n), let row = ordersOfRootsOfUnityInteger !! (n-1) in map (row!!) $ map pred $ take 10 $ iterate (2*) 1)) [2..30] --} - -ordersOfRootsOfUnityIntegerCondensed :: [[Int]] -ordersOfRootsOfUnityIntegerCondensed = - flip map [1..] $ \modu -> - let maxOrder = maximumOrderOfPrimitiveRootsOfUnity (modu::Integer) - in map (length . rootsOfUnityPower modu) $ --- filter (flip divides maxOrder) $ - [Order 1 .. maxOrder] - -rootsOfUnityPower :: - (PrimitiveRoot a, Eq a) => a -> Order -> [a] -rootsOfUnityPower modu (Order expo) = - filter - (\n -> - PID.coprime n modu - && - RC.representative ((n /: modu) ^ expo) == one) $ - primitiveRootCandidates modu - -{- -Corollary from the Carmichael function properties: -If in Z_n there is a primitive root r of unity of order o, -then for every Z_{m \cdot n} there is also a primitive root of unity -with the same order. - -Collary: -If in Z_n1 you have a primitive root of order o1, -and so on for Z_{n_k} and order ok, -then Z_{lcm(n1,...,nk)} has primitive roots for every of the order o1,...,on. - -Conjecture: -If Z_n has a total number of m primitive roots of unity of order o, -then Z_{k*n} has at least m primitive roots of unity of order o. - -Conjecture: -Precondition: In Z_n there is a primitive root of prime order o. -Claims: -a) There are natural numbers m and k with n = m*(k*o+1) or n = m*o. -b) The smallest such n is of the form k*o+1 with k>1. - -Counterexample for a) and non-prime o: o=15, n=77 -Counterexample for b) and non-prime o: - o=20, n=25; o=27, n=81; o=54, n=81; o=55, n=121 - -Corollary from definition of Carmichael function: -For n>1, Z_{2^{n+2}} has primitive root of unity of order 2^n. --} - -{- | -Given an order find integer residue rings -where roots of unity of this order exist. -The way they are constructed also warrants, -that 'order' is a unit (i.e. invertible) in those rings. - -The list is not exhaustive -but computes suggestions quickly. -The first found modulus seems to be smallest one that exist. -However, the first modulus is not the smallest one -among the ones that only have the wanted primitive root, -but where 'order' is not necessarily a unit. -E.g. - -> ringsWithPrimitiveRootOfUnityAndUnit 840 == 2521 : 3361 : ... - -but the smallest modulus is 1763. - -Most of the numbers are primes. -For these the recursion property of the Carmichael function -immediately yields that there are primitive roots of unity of order 'order'. --} -ringsWithPrimitiveRootOfUnityAndUnit :: Order -> [Integer] -ringsWithPrimitiveRootOfUnityAndUnit order@(Order k) = - filter (flip hasPrimitiveRootOfUnityInteger order) $ - iterate (k+) 1 - - -ringsWithPrimitiveRootsOfUnityAndUnitsNaive :: [Order] -> [Integer] -> [Integer] -ringsWithPrimitiveRootsOfUnityAndUnitsNaive rootOrders units = - filter - (\n -> - all (hasPrimitiveRootOfUnityInteger n) rootOrders && - all (PID.coprime n) units) - [1..] - - -{- -It would be nice to have the Omega monad here -in order to enumerate all rings. --} -ringWithPrimitiveRootsOfUnityAndUnits :: [Order] -> [Integer] -> Integer -ringWithPrimitiveRootsOfUnityAndUnits rootOrders units = - let p = lcmMulti units - in lcmMulti $ - map (head . filter (PID.coprime p) . - ringsWithPrimitiveRootOfUnityAndUnit) $ - rootOrders - -{- -Search for an appriopriate modulus -using the recursive definition of the Carmichael function. -We chose the prime factors of the Carmichael function argument -such that we get at least the prime factors in the function value that we need. - -The modulus constructed this way is usually not the smallest possible -and it also does not provide that 'n' is a unit in the residue ring. -The simpler function 'ringsWithPrimitiveRootOfUnityAndUnit' -will usually produce a smaller modulus. --} -ringWithPrimitiveRootsOfUnity :: Order -> Integer -ringWithPrimitiveRootsOfUnity (Order n) = - case n of - 0 -> 2 - _ -> - product $ map (uncurry ringPower) $ snd $ - mapAccumL - (\factors (e,p) -> - if Map.findWithDefault 0 p factors >= e - then (factors, (0,p)) - else - if p==2 - then - (factors, - case e of - 0 -> (0,2) - 1 -> (1,3) - 2 -> (1,5) - _ -> (e+2, 2)) - else - (Map.unionWith max factors $ - Map.fromList $ map swap $ primeFactors $ p-1, - (e+1, p))) - Map.empty $ - reverse $ primeFactors $ lcmMulti $ - n : map (subtract 1) (partialPrimes n) - -lcmMulti :: (PID.C a) => [a] -> a -lcmMulti = foldl lcm one - - -{- | -List all numbers that only contain prime factors 2 and 3 in ascending order. -OEIS:A003586 --} -numbers3Smooth :: [Integer] -numbers3Smooth = - foldr - (\(x0:x1:xs) ys -> x0 : x1 : ListHT.mergeBy (<=) xs ys) - (error "numbers3Smooth: infinite list should not have an end") $ - iterate (map (3*)) $ - iterate (2*) 1 - -numbers3SmoothAlt :: [Integer] -numbers3SmoothAlt = - unfoldr - (fmap (\(m,rest) -> (m, Set.union rest $ Set.fromAscList [2*m,3*m])) . - Set.minView) $ - Set.singleton 1 - -{- -OEIS:A051037 --} -numbers5Smooth :: [Integer] -numbers5Smooth = - foldr - (\(x0:x1:x2:xs) ys -> x0 : x1 : x2 : ListHT.mergeBy (<=) xs ys) - (error "numbers5Smooth: infinite list should not have an end") $ - iterate (map (5*)) $ - numbers3Smooth - -numbers5SmoothAlt :: [Integer] -numbers5SmoothAlt = - unfoldr - (fmap (\(m,rest) -> (m, Set.union rest $ Set.fromAscList [2*m,3*m,5*m])) . - Set.minView) $ - Set.singleton 1 - -ceilingPowerOfTwo :: (Ring.C a, Bits a) => a -> a -ceilingPowerOfTwo 0 = 1 -ceilingPowerOfTwo n = - (1+) $ fst $ head $ - dropWhile (uncurry (/=)) $ - ListHT.mapAdjacent (,) $ - scanl (\m d -> shiftR m d .|. m) (n-1) $ - iterate (2*) 1 - -divideByMaximumPower :: - (Integral.C a, ZeroTestable.C a) => a -> a -> a -divideByMaximumPower b n = - last $ - n : unfoldr (\m -> case divMod m b of (q,r) -> toMaybe (isZero r) (q,q)) n - -divideByMaximumPowerRecursive :: - (Integral.C a, Eq a, ZeroTestable.C a) => a -> a -> a -divideByMaximumPowerRecursive b = - let recourse n = - case divMod b n of - (q,0) -> recourse q - _ -> n - in recourse - -getMaximumExponent :: - (Integral.C a, ZeroTestable.C a) => - a -> a -> (Int,a) -getMaximumExponent b n = - last $ (0,n) : - unfoldr - (\(e,m) -> - let (q,r) = divMod m b - eq = (e+1,q) - in toMaybe (isZero r) (eq,eq)) - (0,n) - -is3Smooth :: Integer -> Bool -is3Smooth = - (1==) . - divideByMaximumPower 3 . - divideByMaximumPower 2 - -is5Smooth :: Integer -> Bool -is5Smooth = - (1==) . - divideByMaximumPower 5 . - divideByMaximumPower 3 . - divideByMaximumPower 2 - -{- | -Compute the smallest composite of 2 and 3 that is as least as large as the input. -This can be interpreted as solving an integer linear programming problem with -min (\(a,b) -> a * log 2 + b * log 3) -over the domain {(a,b) : a>=0, b>=0, a * log 2 + b * log 3 >= log n} --} -{- -Problem: We cannot just start with the ceilingPowerOfTwo -and then multiply with 3/4 until we fall below n, -since the 3/4 decreases too fast. -27/32 is closer to one, -and higher powers of 3 and 2 in the ratio make the ratio even closer to one. --} -ceiling3Smooth :: Integer -> Integer -ceiling3Smooth n = - head $ dropWhile ( Integer -ceiling5Smooth n = - head $ dropWhile ( Integer -ceiling3SmoothNaive = - head . - dropWhile (not . is3Smooth) . - iterate (1+) - -ceiling5SmoothNaive :: Integer -> Integer -ceiling5SmoothNaive = - head . - dropWhile (not . is5Smooth) . - iterate (1+) - - -{- | -Compute all primes that occur in the course of dividing -a Fourier transform into sub-transforms. --} -partialPrimes :: Integer -> [Integer] -partialPrimes = - let primeFactorSet = - Set.fromAscList . uniquePrimeFactors - in unfoldr - (fmap - (\(p,set) -> - (p, Set.union (primeFactorSet (p-1)) set)) . - Set.maxView) - . - primeFactorSet - --- cf. htam:NumberTheory -uniquePrimeFactors :: - (Integral.C a, Bits a, ZeroTestable.C a, Ord a) => - a -> [a] --- map snd . primeFactors -uniquePrimeFactors n = - let oddFactors = - foldr - (\p go m -> - let (q,r) = divMod m p - in if r==0 - then p : go (divideByMaximumPower p q) - else - if q >= p - then go m - else if m==1 then [] else m : []) - (error "uniquePrimeFactors: end of infinite list") - (iterate (2+) 3) - in case powerOfTwoFactors n of - (1,m) -> oddFactors m - (_,m) -> 2 : oddFactors m - -{- | -Prime factors and their exponents in ascending order. --} -primeFactors :: - (PrimitiveRoot a, Ord a) => a -> [(Int, a)] -primeFactors n = - let oddFactors = - foldr - (\p go m -> - let (q0,r) = divMod m p - in if r==0 - then - case getMaximumExponent p q0 of - (e,q1) -> (e+1,p) : go q1 - else - if q0 >= p - then go m - else if m==1 then [] else (1,m) : []) - (const []) - (filter (not . Units.isUnit) $ - primitiveRootCandidates n) - in case getMaximumExponent 2 n of - (0,m) -> oddFactors m - (e,m) -> (e,2) : oddFactors m - -{- -cf. htam:NumberTheory - -Shall this be moved to NumericPrelude? - -It should be replaced by a more sophisticated prime test. --} -isPrime :: Integer -> Bool -isPrime n = - case primeFactors n of - [] -> False - (e,m):_ -> e==1 && m==n - -{- | -Find lengths of signals that require many interim Rader transforms -and end with the given length. - -raderWorstCases 2 = OEIS-A061092 -raderWorstCases 5 = tail OEIS-A059411 - -Smallest raderWorstCase numbers are 2,5,13,17,19,31,37,41,43,61,... -This matches the definition of OEIS-A061303. --} -raderWorstCases :: Integer -> [Integer] -raderWorstCases = - iterate - (\n -> - head $ dropWhile (not . isPrime) $ - tail $ iterate (n+) 1) - -{- | -This is usually faster than 'fastFourierRing' -since it does not need to factor large numbers. -However, the generated modulus is usually much greater. --} -{- -I see the following opportunities for optimization: - -1. Speedup 'fastFourierRing' by - faster primality test (Miller-Rabin) and - faster prime factorization (Pollard-Rho-method). - These are also needed for - maximumOrderOfPrimitiveRootsOfUnityInteger - that is used by Fourier.Element.primitiveRoot - in order to compute a root with maximum order. - -2. Reduce the moduli produced by 'fastFourierRingAlt' - by merging some orders that are passed to - ringWithPrimitiveRootsOfUnityAndUnits, - such that an LCM of a group of orders can still be handled. - This is a kind of knapsack problem. - Maybe we could collect the factors in a way - such that (lcm orderGroup + 1) is prime. - -3. Avoid to compute factorizations of numbers - where we already know the factors - or where we do not need the factors at all. - Use the factors returned by partialPrimes - in order to compute a prime factorization - of lcmMulti (map pred (partialPrimes n)). - Call this (product ps). - Now search for rings of moduli (1 + k * product ps), - where there are (small) primitive roots of order (product ps). - We only need to check whether there are small numbers - such as 2, 3, 5, 6, 7 that have a (product ps)-th power that is 1, - using fast exponentiation. - If there is a power being 1 then check for primitivity - by computing (k * product ps / p)-th powers - for all prime factors p of (k * product ps). - If there is no primitive root <= 7, - there may still be a primitive root of wanted order, - but it is then cheaper to seek for larger moduli. - - If we finally have a nice modulus - it is still stupid to factorize (modulus-1) - and search for a primitive root - in each invocation of Fourier.Element.primitiveRoot. - We could define a special datatype analogously to ResidueClass, - that stores the primitive root and its order - additional to the ResidueClass modulus. --} -fastFourierRingAlt :: Int -> Integer -fastFourierRingAlt n = - case n of - 0 -> 2 - 1 -> 2 - _ -> - let ni = fromIntegral n - ps = filter (>1) (map (subtract 1) (partialPrimes ni)) - in ringWithPrimitiveRootsOfUnityAndUnits (map Order $ ni : ps) ps - -{- | -Determine an integer residue ring -in which a Fast Fourier transform of size n can be performed. -It must contain certain primitive roots. -If we choose a non-primitive root, -then some off-diagonal values in F^-1·F are non-zero. --} -{- -When we need roots of orders o1,...,ok and according inverse elements -we can also ask for a ring, where there is a root of order lcm(o1,...,ok). -The answer to both questions is the same set of rings. -This can be proven using the statement, -that the order of any primitive root -divides the carmichael value of the modulus. - -Since ringWithPrimitiveRootsOfUnityAndUnits -multiplies the moduli of rings for o1,...,ok, -it will produce large moduli. --} -fastFourierRing :: Int -> Integer -fastFourierRing n = - case n of - 0 -> 2 - 1 -> 2 - _ -> - let ni = fromIntegral n - in {- - We cannot use ringsWithPrimitiveRootOfUnityAndUnit - since for 359 we already get an Int overflow. - For 719, 1439, 2879 we also get a very large value. - -} - head $ filter isPrime $ - (\order -> iterate (order +) 1) $ - lcmMulti $ - ni : map (subtract 1) (partialPrimes ni) diff --git a/src/Synthesizer/Generic/Permutation.hs b/src/Synthesizer/Generic/Permutation.hs deleted file mode 100644 --- a/src/Synthesizer/Generic/Permutation.hs +++ /dev/null @@ -1,151 +0,0 @@ -{- | -Permutations of signals as needed for Fast Fourier transforms. -Most functions are independent of the Signal framework. -We could move them as well to Synthesizer.Basic. --} -module Synthesizer.Generic.Permutation where - -import qualified Synthesizer.Basic.NumberTheory as NumberTheory - -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.State.Signal as SigS - -import qualified Data.StorableVector.ST.Strict as SVST -import qualified Data.StorableVector as SV - -import qualified Algebra.PrincipalIdealDomain as PID - - - -type T = SV.Vector Int - -apply :: - (SigG.Transform sig y) => - T -> sig y -> sig y -apply p xs = - SigG.takeStateMatch xs $ - SigS.map (SigG.index xs) $ - SigS.fromStrictStorableSignal p - - -size :: T -> Int -size = SV.length - - -{- | -> inverse (transposition n m) = transposition m n --} -transposition :: - Int -> Int -> T -transposition n m = - fst $ SV.unfoldrN (n*m) - (\(i,j,k0) -> Just (i, - case pred k0 of - 0 -> let j1 = j+1 in (j1, j1, m) - k1 -> (i+n, j, k1))) - (0,0,m) - - -{- -In general the inverse of a skewGrid -does not look like even a generalized skewGrid. -E.g. @inverse $ skewGrid 3 4@. --} -skewGrid :: - Int -> Int -> T -skewGrid n m = - let len = n*m - in fst $ SV.unfoldrN len - (\(i0,k0) -> Just (i0, - let k1 = pred k0 - i1 = i0+n - in if k1==0 - then (mod (i1+m) len, m) - else (mod i1 len, k1))) - (0,m) - -{- | -> inverse (skewGrid n m) == skewGridInv n m - -In general the inverse of a skewGrid -cannot be expressed like skewGrid or skewGridCRT. -E.g. @inverse $ skewGrid 3 4@. --} -skewGridInv :: - Int -> Int -> T -skewGridInv n m = - SV.pack $ - map (\k -> - let Just (i,j) = PID.diophantine k n m - in mod i m + mod j n * m) $ - take (n*m) $ iterate (1+) 0 - -skewGridCRT :: - Int -> Int -> T -skewGridCRT n m = - let len = n*m - (ni,mi) = snd $ PID.extendedGCD n m - in fst $ SV.unfoldrN len - (\(i0,k0) -> Just (i0, - let k1 = pred k0 - i1 = i0+ni*n - in if k1==0 - then (mod (i1+mi*m) len, m) - else (mod i1 len, k1))) - (0,m) - -skewGridCRTInv :: - Int -> Int -> T -skewGridCRTInv n m = - fst $ SV.packN (n*m) $ - map (\k -> mod k m + mod k n * m) $ - iterate (1+) 0 - - -{- | -Beware of 0-based indices stored in the result vector. --} -multiplicative :: Int -> T -multiplicative ni = - let n = fromIntegral ni - gen = NumberTheory.multiplicativeGenerator n - in {- - Since 'gen' is usually 2 or 3, - the error should occur really only for huge signals. - -} - if gen * n > fromIntegral (maxBound :: Int) - then error "signal too long for Int indexing" - else fst $ SV.unfoldrN (ni-1) - (\x -> Just (x-1, mod (fromInteger gen * x) ni)) 1 - -{- | -We only need to compute the inverse permutation explicitly, -because not all signal structures support write to arbitrary indices, -thus Generic.Write does not support it. -For strict StorableVector it would be more efficient -to build the vector directly. - -It holds: - -> inverse . inverse == id --} -inverse :: T -> T -inverse perm = - SVST.runSTVector - (do inv <- SVST.new_ (SV.length perm) - SigS.sequence_ $ - SigS.zipWith (SVST.write inv) - (SigS.fromStrictStorableSignal perm) - (SigS.iterate (1+) 0) - return inv) - -reverse :: T -> T -reverse perm = - fst $ SV.unfoldrN (SV.length perm) - (\mn -> Just $ - case mn of - Nothing -> (SV.head perm, Just $ SV.length perm) - Just n -> - let n1 = n-1 - in (SV.index perm n1, Just n1)) - Nothing diff --git a/src/Test/Main.hs b/src/Test/Main.hs deleted file mode 100644 --- a/src/Test/Main.hs +++ /dev/null @@ -1,49 +0,0 @@ -module Main where - -import qualified Test.Sound.Synthesizer.Plain.Analysis as Analysis -import qualified Test.Sound.Synthesizer.Plain.Control as Control -import qualified Test.Sound.Synthesizer.Plain.Filter as Filter -import qualified Test.Sound.Synthesizer.Plain.Interpolation as Interpolation -import qualified Test.Sound.Synthesizer.Plain.Oscillator as Oscillator -import qualified Test.Sound.Synthesizer.Plain.Wave as Wave -import qualified Test.Sound.Synthesizer.Basic.NumberTheory as NumberTheory -import qualified Test.Sound.Synthesizer.Basic.ToneModulation as ToneModulation -import qualified Test.Sound.Synthesizer.Plain.ToneModulation as ToneModulationL -import qualified Test.Sound.Synthesizer.Generic.ToneModulation as ToneModulationG -import qualified Test.Sound.Synthesizer.Generic.Permutation as Permutation -import qualified Test.Sound.Synthesizer.Generic.Fourier as Fourier -import qualified Test.Sound.Synthesizer.Generic.FourierInteger as FourierInteger -import qualified Test.Sound.Synthesizer.Generic.Filter as FilterG -import qualified Test.Sound.Synthesizer.Generic.Cut as CutG -import qualified Test.Sound.Synthesizer.Causal.Analysis as AnalysisC -import qualified Test.Sound.Synthesizer.Storable.Cut as Cut - -import Data.Tuple.HT (mapFst, ) - - -prefix :: String -> [(String, IO ())] -> [(String, IO ())] -prefix msg = - map (mapFst (\str -> msg ++ "." ++ str)) - -main :: IO () -main = - mapM_ (\(msg,io) -> putStr (msg++": ") >> io) $ - concat $ - prefix "Plain.Analysis" Analysis.tests : - prefix "Plain.Control" Control.tests : - prefix "Plain.Filter" Filter.tests : - prefix "Plain.Interpolation" Interpolation.tests : - prefix "Plain.Oscillator" Oscillator.tests : - prefix "Plain.Wave" Wave.tests : - prefix "Storable.Cut" Cut.tests : - prefix "Generic.Cut" CutG.tests : - prefix "Basic.ToneModulation" ToneModulation.tests : - prefix "Plain.ToneModulation" ToneModulationL.tests : - prefix "Generic.ToneModulation" ToneModulationG.tests : - prefix "Generic.Permutation" Permutation.tests : - prefix "Generic.Fourier" Fourier.tests : - prefix "Basic.NumberTheory" NumberTheory.tests : - prefix "Generic.FourierInteger" FourierInteger.tests : - prefix "Generic.Filter" FilterG.tests : - prefix "Causal.Analysis" AnalysisC.tests : - [] diff --git a/src/Test/Sound/Synthesizer/Basic/NumberTheory.hs b/src/Test/Sound/Synthesizer/Basic/NumberTheory.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Basic/NumberTheory.hs +++ /dev/null @@ -1,119 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Sound.Synthesizer.Basic.NumberTheory (tests) where - -import Synthesizer.Basic.NumberTheory (Order(Order), ) -import qualified Synthesizer.Basic.NumberTheory as NT -import qualified Data.Set as Set - -import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) - -import qualified Algebra.Absolute as Absolute - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -newtype Cardinal a = Cardinal a - deriving (Show) - -instance (Absolute.C a, Arbitrary a) => Arbitrary (Cardinal a) where - arbitrary = fmap (Cardinal . abs) arbitrary - - -newtype Positive a = Positive a - deriving (Show) - -instance (Absolute.C a, Arbitrary a) => Arbitrary (Positive a) where - arbitrary = fmap (Positive . (1+) . abs) arbitrary - - -simple :: - (Testable t, Arbitrary (wrapper Integer), Show (wrapper Integer)) => - (wrapper Integer -> t) -> IO () -simple = quickCheck - -tests :: [(String, IO ())] -tests = - ("primitiveRootsOfUnity naive vs. power", - simple $ \(Cardinal m) order -> - NT.primitiveRootsOfUnityNaive m order - == - NT.primitiveRootsOfUnityPower m order) : - ("primitiveRootsOfUnity naive vs. fullorbit", - simple $ \(Cardinal m) order -> - NT.primitiveRootsOfUnityNaive m order - == - (Set.toAscList $ Set.fromList $ - NT.primitiveRootsOfUnityFullOrbit m order)) : - ("Carmichael theorem", - simple $ \(Positive a) (Positive b) -> - NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity (lcm a b)) - == - lcm - (NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity a)) - (NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity b))) : - ("maximumOrderOfPrimitiveRootsOfUnity naive vs. integer", - simple $ \(Positive m) -> - NT.maximumOrderOfPrimitiveRootsOfUnityNaive m - == - NT.maximumOrderOfPrimitiveRootsOfUnityInteger m) : - ("number of rootsOfUnityPower, lcm", - simple $ \(Positive m) ao@(Order a) bo@(Order b) -> - let g = length . NT.rootsOfUnityPower m - in g (Order $ lcm a b) == lcm (g ao) (g bo)) : - ("ringsWithPrimitiveRootsOfUnityAndUnits: minimal modulus", - quickCheck $ \order@(Order expo) -> - (head $ NT.ringsWithPrimitiveRootOfUnityAndUnit order) - == - (head $ NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive - [order] [expo])) : - ("combine two rings with primitive roots of certain orders", - quickCheck $ \m n -> - let r = lcm - (head (NT.ringsWithPrimitiveRootOfUnityAndUnit m)) - (head (NT.ringsWithPrimitiveRootOfUnityAndUnit n)) - in NT.hasPrimitiveRootOfUnityInteger r m - && - NT.hasPrimitiveRootOfUnityInteger r n) : - ("combine many rings with primitive roots of certain orders", - quickCheck $ \n0 ns0 -> - let ns = take 3 $ map (\n -> 1 + mod n 10) (n0:ns0) - order = NT.lcmMulti ns - in take 3 (NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive - (map Order ns) ns) - == - take 3 (NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive - [Order order] [order])) : -{- -Unfortunately rings with certain units cannot be combined -while maintaining these elements as units. - -Counterexample: - ringsWithPrimitiveRootOfUnityAndUnit 2 = 3:... - ringsWithPrimitiveRootOfUnityAndUnit 3 = 7:... - But in Z_{3·7} the number 3 is no unit. - - ("combine rings with certain units", - quickCheck $ \(Positive m) (Positive n) -> - let r = fromIntegral $ lcm - (head (NT.ringsWithPrimitiveRootOfUnityAndUnit m)) - (head (NT.ringsWithPrimitiveRootOfUnityAndUnit n)) - in PID.coprime r m && PID.coprime r n) : --} - ("number of roots of unity lcm", - quickCheck $ \(Positive n) (Positive k) (Positive l) -> - let orders = NT.ordersOfRootsOfUnityInteger !! (n-1) - in lcm (orders!!(k-1)) (orders!!(l-1)) - == - orders !! (lcm k l - 1)) : - ("number of roots of unity vs. primitive roots", - quickCheck $ \(Positive n) (Positive k) -> - (sum $ map snd $ - filter (flip divides k . fst) $ - zip - [1..] - (NT.ordersOfPrimitiveRootsOfUnityInteger !! (n-1))) - == - NT.ordersOfRootsOfUnityInteger !! (n-1) !! (k-1)) : - [] diff --git a/src/Test/Sound/Synthesizer/Basic/ToneModulation.hs b/src/Test/Sound/Synthesizer/Basic/ToneModulation.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Basic/ToneModulation.hs +++ /dev/null @@ -1,93 +0,0 @@ -module Test.Sound.Synthesizer.Basic.ToneModulation where - -import qualified Synthesizer.Interpolation as Interpolation -import Synthesizer.Interpolation (margin, ) - -import qualified Synthesizer.Basic.Phase as Phase -import qualified Synthesizer.Basic.ToneModulation as ToneMod - -import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest - -import Test.QuickCheck (quickCheck, Property, (==>), Testable, ) --- import Test.Utility - -import qualified Number.NonNegative as NonNeg - -import qualified Algebra.RealField as RealField -import qualified Algebra.Field as Field - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -untangleShapePhase :: (Field.C a, Eq a) => - Int -> a -> (a, a) -> Property -untangleShapePhase periodInt period c = - period /= zero ==> - ToneMod.untangleShapePhase periodInt period c == - ToneMod.untangleShapePhaseAnalytic periodInt period c - -flattenShapePhase :: (RealField.C a) => - Int -> a -> (a, Phase.T a) -> Property -flattenShapePhase periodInt period c = - period /= zero ==> - ToneMod.flattenShapePhase periodInt period c == - ToneMod.flattenShapePhaseAnalytic periodInt period c - - --- * auxiliary quickCheck functions - -{- -Although that looks like a too small value, it is actually right, -because numberLeap counts intervals of size periodInt, not single elements. -So numberLeap=2 like in linear interpolation means 2*periodInt. --} -minLength :: - Interpolation.T a v -> - Interpolation.T a v -> - Int -> NonNeg.Int -> Int -minLength ipLeap ipStep = - minLengthMargin (margin ipLeap) (margin ipStep) - -minLengthMargin :: - Interpolation.Margin -> - Interpolation.Margin -> - Int -> NonNeg.Int -> Int -minLengthMargin marginLeap marginStep periodInt ext = - ToneMod.interpolationNumber - marginLeap marginStep periodInt + - NonNeg.toNumber ext - - - -shapeLimits :: - Interpolation.T a v -> - Interpolation.T a v -> - Int -> Int -> (Int, Int) -shapeLimits ipLeap ipStep periodInt len = - ToneMod.shapeLimits - (margin ipLeap) (margin ipStep) - periodInt len - - - -testRationalLineIp :: Testable quickCheck => - (InterpolationTest.LinePreserving Rational Rational -> quickCheck) -> IO () -testRationalLineIp f = quickCheck f - -testRationalIp :: Testable quickCheck => - (InterpolationTest.T Rational Rational -> quickCheck) -> IO () -testRationalIp f = quickCheck f - - -tests :: [(String, IO ())] -tests = - ("untangleShapePhase", - quickCheck $ \periodInt period -> - untangleShapePhase periodInt (period :: Rational)) : - ("flattenShapePhase", - quickCheck $ \periodInt period -> - flattenShapePhase periodInt (period :: Rational)) : - [] diff --git a/src/Test/Sound/Synthesizer/Causal/Analysis.hs b/src/Test/Sound/Synthesizer/Causal/Analysis.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Causal/Analysis.hs +++ /dev/null @@ -1,32 +0,0 @@ -module Test.Sound.Synthesizer.Causal.Analysis (tests) where - -import qualified Synthesizer.Causal.Analysis as AnaC -import qualified Synthesizer.Causal.Process as Causal -import qualified Synthesizer.Plain.Analysis as Ana - -import Control.Arrow ((<<<), ) - -import qualified Data.List.Match as Match - -import Test.QuickCheck (quickCheck, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -tests :: [(String, IO ())] -tests = - ("deltaSigmaModulation", - quickCheck $ \xs -> - Match.take xs (Ana.deltaSigmaModulation xs) - == - Causal.apply AnaC.deltaSigmaModulation (xs::[Rational])) : - ("deltaSigmaModulationPositive", - quickCheck $ \threshold xs -> - Match.take xs (Ana.deltaSigmaModulationPositive threshold xs) - == - Causal.apply - (AnaC.deltaSigmaModulationPositive <<< - Causal.feedConstFst threshold) (xs::[Rational])) : - [] diff --git a/src/Test/Sound/Synthesizer/Generic/Cut.hs b/src/Test/Sound/Synthesizer/Generic/Cut.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/Cut.hs +++ /dev/null @@ -1,104 +0,0 @@ -module Test.Sound.Synthesizer.Generic.Cut (tests) where - -import qualified Synthesizer.Generic.Cut as CutG -import qualified Synthesizer.Generic.Signal as SigG - -import qualified Synthesizer.Storable.Signal as SigSt - -import qualified Synthesizer.ChunkySize.Signal as SigChunky -import qualified Synthesizer.ChunkySize as ChunkySize - -import qualified Data.StorableVector as SV -import qualified Data.StorableVector.Lazy.Pattern as SVP - -import qualified Synthesizer.State.Signal as SigS - -import qualified Data.EventList.Relative.BodyTime as EventList - -import qualified Number.NonNegative as NonNeg -import qualified Number.NonNegativeChunky as Chunky - -import qualified Numeric.NonNegative.Wrapper as NonNeg98 - -import Data.Tuple.HT (mapSnd, ) - -import Test.QuickCheck (quickCheck, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -dropMarginRemLength :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool -dropMarginRemLength nn nm xs = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - (k,ys) = CutG.dropMarginRem n m xs - in length xs - m == length ys - k - -dropMarginRemState :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool -dropMarginRemState nn nm xs = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - in CutG.dropMarginRem n m (SigS.fromList xs) - == - mapSnd SigS.fromList (CutG.dropMarginRem n m xs) - -dropMarginRemSV :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool -dropMarginRemSV nn nm xs = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - in CutG.dropMarginRem n m (SV.pack xs) - == - mapSnd SV.pack (CutG.dropMarginRem n m xs) - -dropMarginRemSVL :: NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> [Int] -> Bool -dropMarginRemSVL nn nm pat xs = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - in CutG.dropMarginRem n m - (CutG.take (CutG.length pat) xs) - == - mapSnd SigG.toList - (CutG.dropMarginRem n m - (SigChunky.fromState pat $ - SigG.toState xs :: SigSt.T Int)) - -dropMarginRemChunkySize :: - NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> Int -> Bool -dropMarginRemChunkySize nn nm pat x = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - in CutG.dropMarginRem n m pat - == - mapSnd - (ChunkySize.fromStorableVectorSize . SVP.length) - (CutG.dropMarginRem n m - (SVP.replicate (ChunkySize.toStorableVectorSize pat) x)) - -dropMarginRemPiecewise :: - NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> Int -> Bool -dropMarginRemPiecewise nn nm pat x = - let n = NonNeg.toNumber nn - m = NonNeg.toNumber nm - in CutG.dropMarginRem n m pat - == - mapSnd - (Chunky.fromChunks . - map (\size -> SigG.LazySize $ NonNeg98.toNumber size) . - EventList.getTimes) - (CutG.dropMarginRem n m - (EventList.fromPairList $ map ((,) x) $ - map (\(SigG.LazySize size) -> NonNeg98.fromNumber size) $ - Chunky.toChunks pat)) - - -tests :: [(String, IO ())] -tests = - ("dropMarginRemLength", quickCheck dropMarginRemLength) : - ("dropMarginRemState", quickCheck dropMarginRemState) : - ("dropMarginRemSV", quickCheck dropMarginRemSV) : - ("dropMarginRemSVL", quickCheck dropMarginRemSVL) : - ("dropMarginRemChunkySize", quickCheck dropMarginRemChunkySize) : - ("dropMarginRemPiecewise", quickCheck dropMarginRemPiecewise) : - [] diff --git a/src/Test/Sound/Synthesizer/Generic/Filter.hs b/src/Test/Sound/Synthesizer/Generic/Filter.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/Filter.hs +++ /dev/null @@ -1,64 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Sound.Synthesizer.Generic.Filter (tests) where - -import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG -import qualified Synthesizer.Generic.Cyclic as Cyclic -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.Generic.Cut as CutG -import qualified Synthesizer.Plain.Signal as Sig - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty - -import Test.QuickCheck (Testable, quickCheck, ) - --- import qualified Algebra.Ring as Ring - -import qualified Algebra.Laws as Law - -import NumericPrelude.Numeric -import NumericPrelude.Base - - -simple :: - (Testable t) => - (Sig.T Int -> t) -> IO () -simple = quickCheck - -(=|=) :: - (Eq sig, CutG.Transform sig) => - sig -> sig -> Bool -x =|= y = - CutG.take 100 x == CutG.take 100 y - -tests :: [(String, IO ())] -tests = - ("identity", - simple $ Law.identity FiltNRG.generic $ SigG.singleton one) : - ("commutativity", - simple $ Law.commutative FiltNRG.generic) : - ("distributivity", - simple $ Law.leftDistributive FiltNRG.generic SigG.mix) : - ("karatsuba finite", - simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaFinite (*) x y) : - ("karatsuba finite-infinite", - simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaFiniteInfinite (*) x y) : - ("karatsuba infinite", - simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaInfinite (*) x y) : - ("karatsuba finite-infinite cycle", - simple $ \x yn -> - case NonEmpty.toInfiniteList yn of - y -> FiltNRG.generic x y =|= FiltNRG.karatsubaFiniteInfinite (*) x y) : - ("karatsuba infinite cycle", - simple $ \x yn -> - case NonEmpty.toInfiniteList yn of - y -> FiltNRG.generic x y =|= FiltNRG.karatsubaInfinite (*) x y) : - ("convolve triple", - quickCheck $ \x y -> - Cyclic.sumAndConvolveTriple x y == - Cyclic.sumAndConvolveTripleAlt x (y :: Cyclic.Triple Integer)) : - ("periodic summation", - simple $ \x y n -> - let periodic = Cyclic.fromSignal SigG.defaultLazySize (1 + abs n) - in Cyclic.convolve (periodic x) (periodic y) == - periodic (FiltNRG.generic x y)) : - [] diff --git a/src/Test/Sound/Synthesizer/Generic/Fourier.hs b/src/Test/Sound/Synthesizer/Generic/Fourier.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/Fourier.hs +++ /dev/null @@ -1,151 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Sound.Synthesizer.Generic.Fourier (tests) where - -import qualified Synthesizer.Generic.Fourier as Fourier -import qualified Synthesizer.Generic.Cyclic as Cyclic -import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG -import qualified Synthesizer.Generic.Analysis as AnaG -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.Generic.Cut as CutG -import qualified Synthesizer.Storable.Signal as SigSt -import qualified Synthesizer.State.Signal as SigS - -import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) -import Test.Utility (approxEqualAbs, approxEqualComplexAbs, ) - -import qualified Number.Complex as Complex - -import qualified Algebra.Ring as Ring -import qualified Algebra.Additive as Additive - -import Control.Monad (liftM2, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -tolerance :: Double -tolerance = 1e-10 - -normalize :: - SigSt.T (Complex.T Double) -> SigSt.T (Complex.T Double) -normalize xs = - FiltNRG.amplifyVector - (recip $ max (0.1::Double) $ AnaG.volumeVectorMaximum xs) xs - -newtype Normed = Normed (SigSt.T (Complex.T Double)) - deriving (Show) - -instance Arbitrary Normed where - arbitrary = fmap (Normed . normalize) arbitrary - - -data Normed2 = - Normed2 - (SigSt.T (Complex.T Double)) - (SigSt.T (Complex.T Double)) - deriving (Show) - -instance Arbitrary Normed2 where - arbitrary = - liftM2 - (\x y -> - let len = min (CutG.length x) (CutG.length y) - in Normed2 - (normalize $ CutG.take len x) - (normalize $ CutG.take len y)) - arbitrary - arbitrary - - --- could be moved to NumericPrelude -class Complex a where - conjugate :: a -> a - -instance (Additive.C a) => Complex (Complex.T a) where - conjugate = Complex.conjugate - -scalarProduct :: - (SigG.Read sig y, Ring.C y, Complex y) => - sig y -> sig y -> y -scalarProduct xs ys = - SigS.sum $ - SigS.zipWith (*) - (SigG.toState xs) - (SigS.map conjugate $ SigG.toState ys) - -(=~=) :: - SigSt.T (Complex.T Double) -> - SigSt.T (Complex.T Double) -> - Bool -(=~=) xs ys = - SigG.length xs == SigG.length ys && - (SigG.foldR (&&) True $ - SigG.zipWith (approxEqualComplexAbs tolerance) xs ys) - -simple :: - (Testable t) => - (SigSt.T (Complex.T Double) -> t) -> IO () -simple = quickCheck - -tests :: [(String, IO ())] -tests = - ("fourier inverse", - quickCheck $ \(Normed x) -> - x =~= - (FiltNRG.amplify (recip $ fromIntegral $ SigG.length x) $ - Fourier.transformBackward $ Fourier.transformForward x)) : - ("double fourier = reverse", - quickCheck $ \(Normed x) -> - x =~= - (Cyclic.reverse $ - FiltNRG.amplify (recip $ fromIntegral $ SigG.length x) $ - Fourier.transformForward $ - Fourier.transformForward x)) : - ("fourier of reverse", - quickCheck $ \(Normed x) -> - Cyclic.reverse (Fourier.transformForward x) =~= - Fourier.transformForward (Cyclic.reverse x)) : - ("fourier of conjugate", - quickCheck $ \(Normed x) -> - (SigG.map Complex.conjugate $ Fourier.transformForward x) - =~= - (Fourier.transformForward $ - SigG.map Complex.conjugate $ Cyclic.reverse x)) : - ("additivity", - quickCheck $ \(Normed2 x y) -> - SigG.mix (Fourier.transformForward x) (Fourier.transformForward y) - =~= - Fourier.transformForward (SigG.mix x y)) : - ("isometry", - simple $ \xs x0 -> - let x = normalize (SigG.cons x0 xs) - in approxEqualAbs tolerance - (AnaG.volumeVectorEuclideanSqr $ Fourier.transformForward x) - (fromIntegral (SigG.length x) * - AnaG.volumeVectorEuclideanSqr x)) : - ("unitarity", - quickCheck $ \(Normed2 x y) -> - approxEqualComplexAbs tolerance - (scalarProduct - (Fourier.transformForward x) (Fourier.transformForward y)) - (fromIntegral (SigG.length x) * scalarProduct x y)) : - ("convolution", - quickCheck $ \(Normed2 x y) -> - SigG.zipWith (*) - (Fourier.transformForward x) - (Fourier.transformForward y) - =~= - Fourier.transformForward (Cyclic.convolve x y)) : - ("convolution cyclic", - quickCheck $ \(Normed2 x y) -> - Fourier.convolveCyclic x y - =~= - Cyclic.convolve x y) : - ("convolution long", - quickCheck $ \(Normed x) (Normed y) -> - FiltNRG.karatsubaFinite (*) x y - =~= - Fourier.convolveWithWindow (Fourier.window x) y) : - [] diff --git a/src/Test/Sound/Synthesizer/Generic/FourierInteger.hs b/src/Test/Sound/Synthesizer/Generic/FourierInteger.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/FourierInteger.hs +++ /dev/null @@ -1,178 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Sound.Synthesizer.Generic.FourierInteger (tests) where - -import qualified Synthesizer.Generic.Fourier as Fourier -import qualified Synthesizer.Generic.Cyclic as Cyclic -import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.Generic.Cut as CutG -import qualified Synthesizer.State.Signal as SigS -import qualified Synthesizer.Plain.Signal as Sig - -import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) - -import qualified Synthesizer.Basic.NumberTheory as NT - -import qualified Number.ResidueClass.Check as RC -import Number.ResidueClass.Check ((/:), ) - -import qualified Algebra.ToInteger as ToInteger -import qualified Algebra.IntegralDomain as Integral -import qualified Algebra.Ring as Ring - -import Control.Monad (liftM2, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -newtype Modulus a = Modulus a - deriving (Show) - -instance Ring.C a => Arbitrary (Modulus a) where - arbitrary = fmap (Modulus . (2+) . fromInteger) arbitrary - - -data ModularSignal = - ModularSignal (Modulus Integer) (Sig.T (RC.T Integer)) - deriving (Show) - -instance Arbitrary ModularSignal where - arbitrary = - fmap (uncurry ModularSignal . signal) arbitrary - - -data ModularSignal2 = - ModularSignal2 - (Modulus Integer) (Sig.T (RC.T Integer)) (Sig.T (RC.T Integer)) - deriving (Show) - -instance Arbitrary ModularSignal2 where - arbitrary = - liftM2 - (\x y -> - let len = min (CutG.length x) (CutG.length y) - m = NT.fastFourierRing len - in ModularSignal2 - (Modulus m) - (fmap (/: m) $ CutG.take len x) - (fmap (/: m) $ CutG.take len y)) - arbitrary - arbitrary - -scalarProduct :: - Modulus Integer -> - Sig.T (RC.T Integer) -> Sig.T (RC.T Integer) -> - RC.T Integer -scalarProduct (Modulus m) xs ys = - SigS.foldL (+) (RC.zero m) $ - SigS.zipWith (*) - (SigG.toState xs) - (SigG.toState ys) - -{- -signal :: - Integral.C a => - Modulus a -> Sig.T a -> Sig.T (RC.T a) -signal (Modulus a) = fmap (/: a) --} - -signal :: - Sig.T Integer -> (Modulus Integer, Sig.T (RC.T Integer)) -signal xs = - let m = NT.fastFourierRing $ length xs - in (Modulus m, fmap (/: m) xs) - -modular :: - (Integral.C a, ToInteger.C b) => - Modulus a -> b -> RC.T a -modular (Modulus m) = - RC.fromRepresentative m . fromIntegral - - -simple :: - (Testable t) => - (Sig.T Integer -> t) -> IO () -simple = quickCheck - -tests :: [(String, IO ())] -tests = - ("fourier inverse", - quickCheck $ \(ModularSignal m x) -> - (Fourier.transformBackward $ Fourier.transformForward x) - == - FiltNRG.amplify (modular m $ length x) x) : - ("double fourier = reverse", - quickCheck $ \(ModularSignal m x) -> - (Cyclic.reverse $ - Fourier.transformForward $ - Fourier.transformForward x) - == - FiltNRG.amplify (modular m $ length x) x) : - ("fourier of reverse", - quickCheck $ \(ModularSignal _m x) -> - Cyclic.reverse (Fourier.transformForward x) == - Fourier.transformForward (Cyclic.reverse x)) : - ("homogenity", - quickCheck $ \(ModularSignal m x) y -> - (FiltNRG.amplify (modular m (y::Integer)) $ - Fourier.transformForward x) - == - (Fourier.transformForward $ - FiltNRG.amplify (modular m y) x)) : - ("additivity", - quickCheck $ \(ModularSignal2 _m x y) -> - SigG.mix (Fourier.transformForward x) (Fourier.transformForward y) - == - Fourier.transformForward (SigG.mix x y)) : -{- - ("isometry", - simple $ \xs x0 -> - let (m,x) = signal (SigG.cons x0 xs) - in (AnaG.volumeVectorEuclideanSqr $ Fourier.transformForward x) - == - (modular m (SigG.length x) * - AnaG.volumeVectorEuclideanSqr x)) : --} - ("unitarity", - quickCheck $ \(ModularSignal2 m x y) -> - {- - since there is no equivalent of a complex conjugate - we have to take the scalar product with the backwards transform. - -} - scalarProduct m - (Fourier.transformForward x) (Fourier.transformBackward y) - == - modular m (length x) * scalarProduct m x y) : - ("convolution", - quickCheck $ \(ModularSignal2 _m x y) -> - SigG.zipWith (*) - (Fourier.transformForward x) - (Fourier.transformForward y) - == - Fourier.transformForward (Cyclic.convolve x y)) : - ("convolution cyclic", - quickCheck $ \(ModularSignal2 _m x y) -> - Fourier.convolveCyclic x y - == - Cyclic.convolve x y) : - ("convolution long", - simple $ \x0 y0 -> - let m = Modulus $ NT.fastFourierRing $ - 2 * (NT.ceilingPowerOfTwo $ length x0) - x = fmap (modular m) x0 - y = fmap (modular m) y0 - in fmap (modular m) (FiltNRG.karatsubaFinite (*) x0 y0) - == - Fourier.convolveWithWindow (Fourier.window x) y) : - ("convolution long modular", - simple $ \x0 y0 -> - let m = Modulus $ NT.fastFourierRing $ - 2 * (NT.ceilingPowerOfTwo $ length x0) - x = fmap (modular m) x0 - y = fmap (modular m) (y0 :: Sig.T Integer) - in FiltNRG.karatsubaFinite (*) x y - == - Fourier.convolveWithWindow (Fourier.window x) y) : - [] diff --git a/src/Test/Sound/Synthesizer/Generic/Permutation.hs b/src/Test/Sound/Synthesizer/Generic/Permutation.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/Permutation.hs +++ /dev/null @@ -1,45 +0,0 @@ -{- -wish list: - - custom Permutation type with Arbitrary instance --} -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Sound.Synthesizer.Generic.Permutation (tests) where - -import qualified Synthesizer.Generic.Permutation as Permutation - -import Test.QuickCheck (quickCheck, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -tests :: [(String, IO ())] -tests = - ("inverse transposition", - quickCheck $ \n0 m0 -> - let n = mod n0 100 - m = mod m0 100 - in Permutation.inverse (Permutation.transposition n m) - == - Permutation.transposition m n) : - ("inverse skewGrid", - quickCheck $ \n0 m0 -> - let g = gcd n0 m0 - (n,m) = if g==0 then (0,0) else (abs (div n0 g), abs (div m0 g)) - in Permutation.inverse (Permutation.skewGrid n m) - == - Permutation.skewGridInv n m) : - ("inverse skewGridCRT", - quickCheck $ \n0 m0 -> - let g = gcd n0 m0 - (n,m) = if g==0 then (0,0) else (abs (div n0 g), abs (div m0 g)) - in Permutation.inverse (Permutation.skewGridCRT n m) - == - Permutation.skewGridCRTInv n m) : - {- - reverse (multiplicative (generator n) n) - == - multiplicative (recip $ generator n) n - -} - [] diff --git a/src/Test/Sound/Synthesizer/Generic/ToneModulation.hs b/src/Test/Sound/Synthesizer/Generic/ToneModulation.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Generic/ToneModulation.hs +++ /dev/null @@ -1,304 +0,0 @@ -module Test.Sound.Synthesizer.Generic.ToneModulation (tests) where - -import Test.Sound.Synthesizer.Basic.ToneModulation ( - minLength, - minLengthMargin, --- shapeLimits, --- testRationalLineIp, - testRationalIp, - ) - -import qualified Synthesizer.Causal.ToneModulation as ToneModC -import qualified Synthesizer.Generic.Wave as WaveG - -import qualified Synthesizer.Plain.Signal as Sig -import qualified Synthesizer.Plain.Oscillator as Osci -import qualified Synthesizer.Plain.Interpolation as Interpolation -import qualified Synthesizer.Plain.ToneModulation as ToneModL -import qualified Synthesizer.Plain.Wave as WaveL -import Synthesizer.Interpolation (marginNumber, ) - -import qualified Synthesizer.Causal.Oscillator as OsciC -import qualified Synthesizer.Causal.Process as Causal - -import qualified Synthesizer.State.Signal as SigS - -import qualified Synthesizer.Basic.Wave as Wave -import qualified Synthesizer.Basic.Phase as Phase - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty -import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest - -import Test.QuickCheck (quickCheck, Property, (==>), ) -import Test.Utility (ArbChar, ) --- import Debug.Trace (trace, ) - -import qualified Number.NonNegative as NonNeg - -import qualified Algebra.RealField as RealField - - -import Data.List.HT (viewL, takeWhileJust, ) -import Data.Tuple.HT (mapSnd, ) -import qualified Data.List as List - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -limitMinRelativeValues :: - Int -> Int -> [NonNeg.Int] -> Bool -limitMinRelativeValues xMin x0 xsnn = - let xs = map NonNeg.toNumber xsnn - (y0,limiter) = ToneModC.limitMinRelativeValues xMin x0 - in (y0, Causal.apply limiter xs) == - ToneModL.limitMinRelativeValues xMin x0 xs - -integrateFractional :: (RealField.C t) => - NonNeg.T t -> t -> Phase.T t -> [NonNeg.T t] -> [t] -> Property -integrateFractional - periodNN shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - (c0, coordinator) = - ToneModC.integrateFractional - period (shape0, phase) - coords = - ToneModL.integrateFractional - period (shape0, shapes) (phase, freqs) - in period /= zero ==> - c0 : Causal.apply coordinator (zip shapes freqs) == - coords - --- oscillatorCellSize :: (Show t, Show v, RealField.C t, Eq v) => -oscillatorCellSize :: (RealField.C t, Eq v) => - Interpolation.Margin -> - Interpolation.Margin -> - NonNeg.Int -> NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> - Property -oscillatorCellSize - marginLeap marginStep periodIntNN periodNN ext - ixs shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - len = minLengthMargin marginLeap marginStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledTone = - ToneModC.oscillatorCells - marginLeap marginStep periodInt period tone - (shape0, Phase.fromRepresentative phase) - `Causal.apply` - zip shapes freqs - in period /= zero && - marginNumber marginLeap > zero && - marginNumber marginStep > zero ==> - all - ((\cell -> - Sig.lengthAtLeast (marginNumber marginLeap) cell && - all (Sig.lengthAtLeast (marginNumber marginStep)) - (take (marginNumber marginLeap) cell)) - . SigS.toList . snd) - resampledTone - -oscillatorSuffixes :: (RealField.C t, Eq v) => - Interpolation.Margin -> - Interpolation.Margin -> - NonNeg.Int -> NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> - Property -oscillatorSuffixes - marginLeap marginStep periodIntNN periodNN ext - ixs shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - len = minLengthMargin marginLeap marginStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledToneA = - init $ - map (\(sp,cell) -> - (sp, takeWhileJust . map (fmap fst . viewL) $ cell)) $ - ToneModL.oscillatorSuffixes - marginLeap marginStep periodInt period tone - (shape0, shapes) (Phase.fromRepresentative phase, freqs) - resampledToneB = - ToneModC.oscillatorSuffixes - marginLeap marginStep periodInt period tone - (shape0, Phase.fromRepresentative phase) - `Causal.apply` - zip shapes freqs - in period /= zero && - periodInt /= zero && - marginNumber marginLeap > zero && - marginNumber marginStep > zero ==> - resampledToneA == resampledToneB - -oscillatorCells :: (RealField.C t, Eq v) => - Interpolation.Margin -> - Interpolation.Margin -> - NonNeg.Int -> NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> - Property -oscillatorCells - marginLeap marginStep periodIntNN periodNN ext - ixs shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - len = minLengthMargin marginLeap marginStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledToneA = - init $ map (mapSnd List.transpose) $ - ToneModL.oscillatorCells - marginLeap marginStep periodInt period tone - (shape0, shapes) (Phase.fromRepresentative phase, freqs) - resampledToneB = - map (mapSnd SigS.toList) $ - ToneModC.oscillatorCells - marginLeap marginStep periodInt period tone - (shape0, Phase.fromRepresentative phase) - `Causal.apply` - zip shapes freqs - in period /= zero && - periodInt /= zero && - marginNumber marginLeap > zero && - marginNumber marginStep > zero ==> - resampledToneA == resampledToneB -{- -Margin {marginNumber = 1, marginOffset = 2} -Margin {marginNumber = 5, marginOffset = 5} -3 % 4 -0 -('\DEL',['~','~','"']) --2 % 5 -2 % 5 -[2 % 1,3 % 4] -[-5 % 2,-1 % 2] --} - -{- | -'WaveL.sampledTone' and 'WaveG.sampledTone' -do not only differ in the signal types they process, -but also in the way they order the signal values. -The cells for 'WaveL.sampledTone' are transposed -with respect to 'WaveG.sampledTone'. --} -sampledTone :: (RealField.C a, Eq v) => - InterpolationTest.T a v -> - InterpolationTest.T a v -> - NonNeg.T a -> NonNeg.Int -> NonEmpty.T v -> - a -> Phase.T a -> Property -sampledTone = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs shape phase -> - let period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - in period /= zero ==> - WaveG.sampledTone ipLeap ipStep period tone shape `Wave.apply` phase == - WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` phase - - - -shapeFreqModFromSampledTone :: (RealField.C t, Eq v) => - InterpolationTest.T t v -> - InterpolationTest.T t v -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> Phase.T t -> [NonNeg.T t] -> [t] -> - Property -shapeFreqModFromSampledTone = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs shape0 phase shapesNN freqs -> - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledToneA = - init $ - Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone - shape0 (Phase.toRepresentative phase) shapes freqs - resampledToneB = - OsciC.shapeFreqModFromSampledTone - ipLeap ipStep period tone shape0 phase - `Causal.apply` - zip shapes freqs - in period /= zero ==> - resampledToneA == resampledToneB - - -{- -We have a problem here with the phase distortion signal, -because frequency and shape modulation control signals -are delayed by one element with respect to the phase distortion. -How can we cope with different lengths of the control signals, -without padding the phase control with zeros? -This one did not work: - phaseDistorts = pd:pds - resampledToneA = - Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone - shape0 (Phase.toRepresentative phase) shapes (init phaseDistorts) freqs --} -shapePhaseFreqModFromSampledTone :: (RealField.C t, Eq v) => - InterpolationTest.T t v -> - InterpolationTest.T t v -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> Phase.T t -> [NonNeg.T t] -> (t,[t]) -> [t] -> - Property -shapePhaseFreqModFromSampledTone = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs shape0 phase shapesNN (pd,pds) freqs -> - let period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - shapes = map NonNeg.toNumber shapesNN - phaseDistorts = pd:pds ++ repeat zero - resampledToneA = - init $ - Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone - shape0 (Phase.toRepresentative phase) shapes phaseDistorts freqs - resampledToneB = - OsciC.shapePhaseFreqModFromSampledTone - ipLeap ipStep period tone shape0 phase - `Causal.apply` - zip3 shapes phaseDistorts freqs - in period /= zero ==> - resampledToneA == resampledToneB - - - -tests :: [(String, IO ())] -tests = - ("limitMinRelativeValues", quickCheck limitMinRelativeValues) : - ("integrateFractional", - quickCheck (\period -> integrateFractional (period :: NonNeg.Rational))) : - ("oscillatorCellSize", - quickCheck (\ml ms periodInt period ext ixs -> - oscillatorCellSize ml ms periodInt (period :: NonNeg.Rational) - ext (ixs :: NonEmpty.T ArbChar))) : - ("oscillatorSuffixes", - quickCheck (\ml ms periodInt period ext ixs -> - oscillatorSuffixes ml ms periodInt (period :: NonNeg.Rational) - ext (ixs :: NonEmpty.T ArbChar))) : - ("oscillatorCells", - quickCheck (\ml ms periodInt period ext ixs -> - oscillatorCells ml ms periodInt (period :: NonNeg.Rational) - ext (ixs :: NonEmpty.T ArbChar))) : - ("sampledTone", - testRationalIp sampledTone) : - ("shapeFreqModFromSampledTone", - testRationalIp shapeFreqModFromSampledTone) : - ("shapePhaseFreqModFromSampledTone", - testRationalIp shapePhaseFreqModFromSampledTone) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Analysis.hs b/src/Test/Sound/Synthesizer/Plain/Analysis.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Analysis.hs +++ /dev/null @@ -1,160 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Analysis (tests) where - -import qualified Synthesizer.Plain.Analysis as Analysis - -import qualified Algebra.Algebraic as Algebraic -import qualified Algebra.RealField as RealField -import qualified Algebra.Field as Field -import qualified Algebra.RealRing as RealRing - -import qualified Algebra.NormedSpace.Maximum as NormedMax -import qualified Algebra.NormedSpace.Euclidean as NormedEuc -import qualified Algebra.NormedSpace.Sum as NormedSum - -import qualified MathObj.LaurentPolynomial as LPoly - -import qualified Data.NonEmpty as NonEmpty -import Data.List (genericLength) - -import Test.QuickCheck (quickCheck, Property, (==>)) -import Test.Utility (approxEqual) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -volumeVectorMaximum :: (NormedMax.C y y, RealRing.C y) => [y] -> Bool -volumeVectorMaximum xs = - Analysis.volumeVectorMaximum xs == Analysis.volumeMaximum xs - -volumeVectorEuclidean :: - (NormedEuc.C y y, Algebraic.C y, Eq y) => - NonEmpty.T [] y -> Bool -volumeVectorEuclidean xs = - let ys = NonEmpty.flatten xs - in Analysis.volumeVectorEuclidean ys == Analysis.volumeEuclidean ys - -volumeVectorEuclideanSqr :: - (NormedEuc.Sqr y y, Field.C y, Eq y) => - NonEmpty.T [] y -> Bool -volumeVectorEuclideanSqr xs = - let ys = NonEmpty.flatten xs - in Analysis.volumeVectorEuclideanSqr ys == Analysis.volumeEuclideanSqr ys - -volumeVectorSum :: - (NormedSum.C y y, RealField.C y) => - NonEmpty.T [] y -> Bool -volumeVectorSum xs = - let ys = NonEmpty.flatten xs - in Analysis.volumeVectorSum ys == Analysis.volumeSum ys - - - -bounds :: Ord a => NonEmpty.T [] a -> Bool -bounds xs = - Analysis.bounds xs == (NonEmpty.minimum xs, NonEmpty.maximum xs) - - -spread :: RealField.C a => (a,a) -> Bool -spread b = - sum (map snd (Analysis.spread b)) == one - - -histogramDiscrete :: NonEmpty.T [] Int -> Bool -histogramDiscrete xs = - Analysis.histogramDiscreteArray xs == - Analysis.histogramDiscreteIntMap xs - -withEmptyHistogram :: - (NonEmpty.T [] y -> (Int, [y])) -> - [y] -> (Int, [y]) -withEmptyHistogram f = - maybe (error "no bounds", []) f . NonEmpty.fetch - -histogramDiscreteLength :: [Int] -> Bool -histogramDiscreteLength xs = - sum (snd (withEmptyHistogram Analysis.histogramDiscreteIntMap xs)) - == - length xs - -histogramDiscreteConcat :: [Int] -> [Int] -> Bool -histogramDiscreteConcat xs ys = - let xHist = withEmptyHistogram Analysis.histogramDiscreteIntMap xs - yHist = withEmptyHistogram Analysis.histogramDiscreteIntMap ys - xyHist0 = - LPoly.add - (uncurry LPoly.Cons xHist) - (uncurry LPoly.Cons yHist) - xyHist1 = - uncurry LPoly.Cons - (withEmptyHistogram Analysis.histogramDiscreteIntMap (xs++ys)) - in if null (LPoly.coeffs xyHist0) - then LPoly.coeffs xyHist0 == LPoly.coeffs xyHist1 - else xyHist0 == xyHist1 - - -histogramLinear :: NonEmpty.T [] Int -> Bool -histogramLinear xs = - let ys = fmap fromIntegral xs :: NonEmpty.T [] Double - in Analysis.histogramLinearArray ys == - Analysis.histogramLinearIntMap ys - - -histogramLinearLength :: NonEmpty.T [] Int -> Bool -histogramLinearLength xs = - let ys = fmap fromIntegral xs :: NonEmpty.T [] Double - in approxEqual 1e-10 - (genericLength $ NonEmpty.tail ys) - (sum (snd (Analysis.histogramLinearIntMap ys))) -{- -With eps = 1e-15 - -Falsifiable, after 83 tests: --20 -[32,-41,11,-25,-17,-27,32,-36,7,-36,38] - -Falsifiable, after 78 tests: -10 -[-35,-28,-28,-24,-4,-29,-14,-29,-20,7,33,-2,-14,-4,7,-40,-5,-12] --} - - - -centroid :: (Field.C a, Eq a) => [a] -> Property -centroid xs = - sum xs /= zero ==> - Analysis.centroid xs == Analysis.centroidAlt xs --- Test.QuickCheck.quickCheck (\xs -> sum xs /= 0 Test.QuickCheck.==> propCentroid (xs::[Rational])) - -histogramDCOffset :: NonEmpty.T (NonEmpty.T []) Int -> Property -histogramDCOffset xs = - let x1 = NonEmpty.flatten xs - x = NonEmpty.flatten x1 - (offset, hist) = Analysis.histogramDiscreteArray x1 - in sum x /= 0 ==> - fromIntegral offset + Analysis.centroid (map fromIntegral hist) == - (Analysis.directCurrentOffset (map fromIntegral x) :: Rational) - - -small :: (Functor f) => f Int -> f Int -small = fmap (flip mod 1000) - - -tests :: [(String, IO ())] -tests = - ("volumeVectorMaximum", quickCheck (volumeVectorMaximum :: [Rational] -> Bool)) : - -- quickCheck may fail due to rounding errors, but so far the computation is exactly the same - ("volumeVectorEuclidean", quickCheck (volumeVectorEuclidean :: NonEmpty.T [] Double -> Bool)) : - ("volumeVectorEuclideanSqr", quickCheck (volumeVectorEuclideanSqr :: NonEmpty.T [] Rational -> Bool)) : - ("volumeVectorSum", quickCheck (volumeVectorSum :: NonEmpty.T [] Rational -> Bool)) : - ("bounds", quickCheck (bounds :: NonEmpty.T [] Rational -> Bool)) : - ("spread", quickCheck (spread :: (Rational,Rational) -> Bool)) : - ("histogramDiscrete", quickCheck (histogramDiscrete . small)) : - ("histogramDiscreteLength", quickCheck (histogramDiscreteLength . small)) : - ("histogramDiscreteConcat", quickCheck (\x y -> histogramDiscreteConcat (small x) (small y))) : - ("histogramLinear", quickCheck (histogramLinear . small)) : - ("histogramLinearLength", quickCheck (histogramLinearLength . small)) : - ("centroid", quickCheck (centroid :: [Rational] -> Property)) : - ("histogramDCOffset", quickCheck (histogramDCOffset . small)) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Control.hs b/src/Test/Sound/Synthesizer/Plain/Control.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Control.hs +++ /dev/null @@ -1,112 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Control (tests) where - -import qualified Synthesizer.Plain.Control as Control - -import Test.QuickCheck (quickCheck, Property, (==>)) -import Test.Utility (equalList, approxEqualListAbs, approxEqualListRel, ) - --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive - -import Data.List (transpose) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -linearRing :: Int -> Int -> Bool -linearRing d y0 = --- Control.linear d y0 == Control.linearMultiscale d y0 - all equalList $ take 100 $ transpose $ - Control.linear d y0 : - Control.linearMultiscale d y0 : - Control.linearStable d y0 : - [] - -{- -*Synthesizer.Plain.Control> propLinearApprox (-2/3) 2 -False - -Need a different definition of approximate equality. --} -linearApprox :: Double -> Double -> Bool -linearApprox d y0 = - all (approxEqualListAbs (1e-10 * max (abs d) (abs y0))) $ - take 100 $ transpose $ - Control.linear d y0 : - Control.linearMean d y0 : - Control.linearMultiscale d y0 : - Control.linearStable d y0 : - [] - -linearExact :: Rational -> Rational -> Bool -linearExact d y0 = - all equalList $ take 100 $ transpose $ - Control.linear d y0 : - Control.linearMean d y0 : - Control.linearMultiscale d y0 : - Control.linearStable d y0 : - [] - -{- -Plain.Control.exponential: Falsifiable, after 88 tests: --8.333333333333326e-2 -3.375 - -Plain.Control.exponential: Falsifiable, after 69 tests: -9.090909090909083e-2 --10.0 - -Plain.Control.exponential: Falsifiable, after 73 tests: --0.125 --1.1428571428571428 - -Plain.Control.exponential2: Falsifiable, after 33 tests: --7.692307692307687e-2 -9.5 --} -exponential :: Double -> Double -> Bool -exponential time y0 = - all (approxEqualListRel (1e-10)) $ take 100 $ transpose $ - Control.exponential time y0 : - Control.exponentialMultiscale time y0 : - Control.exponentialStable time y0 : - [] - -exponential2 :: Double -> Double -> Bool -exponential2 time y0 = - all (approxEqualListRel (1e-10)) $ take 100 $ transpose $ - Control.exponential2 time y0 : - Control.exponential2Multiscale time y0 : - Control.exponential2Stable time y0 : - [] - -cosine :: Double -> Double -> Property -cosine t0 t1 = t0/=t1 ==> - all (approxEqualListAbs (1e-10)) $ - take 100 $ transpose $ - Control.cosine t0 t1 : - Control.cosineMultiscale t0 t1 : - Control.cosineStable t0 t1 : - [] - - -cubic :: (Rational, (Rational, Rational)) -> - (Rational, (Rational, Rational)) -> Property -cubic node0 node1 = fst node0 /= fst node1 ==> - take 100 (Control.cubicHermite node0 node1) == - take 100 (Control.cubicHermiteStable node0 node1) - - - -tests :: [(String, IO ())] -tests = - ("linearRing", quickCheck linearRing) : - ("linearApprox", quickCheck linearApprox) : - ("linearExact", quickCheck linearExact) : - ("exponential", quickCheck exponential) : - ("exponential2", quickCheck exponential2) : - ("cosine", quickCheck cosine) : - ("cubic", quickCheck cubic) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Filter.hs b/src/Test/Sound/Synthesizer/Plain/Filter.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Filter.hs +++ /dev/null @@ -1,199 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Filter (tests) where - -import qualified Synthesizer.Plain.Filter.Recursive.MovingAverage as MA -import qualified Synthesizer.Plain.Filter.NonRecursive as FiltNR -import qualified Synthesizer.Plain.Signal as Sig -import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.Storable.Filter.NonRecursive as FiltNRSt -import qualified Synthesizer.Storable.Signal as SigSt -import qualified Synthesizer.Causal.Filter.NonRecursive as FiltNRC -import qualified Synthesizer.Causal.Process as Causal -import qualified Synthesizer.Frame.Stereo as Stereo - -import qualified Data.StorableVector.Lazy.Pattern as VP - -import Foreign.Storable.Tuple () - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty - -import Test.QuickCheck (quickCheck, {- Property, (==>) -}) -import Test.Utility (equalList, ArbChar, ) - --- import qualified Algebra.Module as Module --- import qualified Algebra.RealField as RealField --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive - -import qualified Number.GaloisField2p32m5 as GF -import qualified Number.NonNegative as NonNeg - -import qualified Numeric.NonNegative.Chunky as Chunky - -import qualified Data.List as List -import Data.Tuple.HT (mapPair, ) - --- import Debug.Trace (trace, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -sums :: NonNeg.Int -> Rational -> Sig.T Rational -> Bool -sums nn x0 xs0 = - let n = min (length xs) (1 + NonNeg.toNumber nn) - xs = x0:xs0 - naive = FiltNR.sums n xs - pyramid = FiltNR.sumsPyramid n xs - rec = drop (n-1) $ MA.sumsStaticInt n xs - in -- this checks only for equal prefixes and can easily go wrong, - -- if one list is empty - and $ zipWith3 (\x y z -> x==y && y==z) naive rec pyramid - -- equalList $ naive : pyramid : rec : [] - -sumRange :: NonNeg.Int -> (NonNeg.Int, NonNeg.Int) -> Sig.T Int -> Bool -sumRange nheight (nl,nr) xs = - let wrap n = mod (NonNeg.toNumber n) (length xs + 1) - height = 1 + NonNeg.toNumber nheight - rng = (wrap nl, wrap nr) - pyr = take height (FiltNR.pyramid xs) - pyrSt = - FiltNRSt.pyramid (+) height - (SigSt.fromList SigSt.defaultChunkSize xs) - in equalList $ - FiltNR.sumRange xs rng : - FiltNR.sumRangeFromPyramid pyr rng : - FiltNR.sumRangeFromPyramidRec pyr rng : - FiltNR.sumRangeFromPyramidFoldr pyr rng : - FiltNRG.sumRangeFromPyramid pyrSt rng : - FiltNRG.sumRangeFromPyramidFoldr pyrSt rng : - FiltNRG.sumRangeFromPyramidReverse pyrSt rng : - [] - -getRange :: (NonNeg.Int, NonNeg.Int) -> NonEmpty.T (NonEmpty.T ArbChar) -> Bool -getRange (nl,nr) pyr0 = - let l = NonNeg.toNumber nl - r = NonNeg.toNumber nr - rng = if l<=r then (l,r) else (r,l) - pyr = map NonEmpty.toInfiniteList $ NonEmpty.toList pyr0 - in equalList $ - FiltNR.getRangeFromPyramid pyr rng : - FiltNRG.consumeRangeFromPyramid (:) [] pyr rng : - [] - -sumsPosModulated :: - NonNeg.Int -> Sig.T (NonNeg.Int,NonNeg.Int) -> NonEmpty.T Int -> Bool -sumsPosModulated nheight nctrl xsc = - let ctrl = map (mapPair (NonNeg.toNumber, NonNeg.toNumber)) nctrl - xs = NonEmpty.toInfiniteList xsc - height = min 10 $ NonNeg.toNumber nheight - in -- trace (show (height, ctrl, xsc)) $ - equalList $ - FiltNR.sumsPosModulated ctrl xs : - FiltNR.sumsPosModulatedPyramid height ctrl xs : - FiltNRG.sumsPosModulatedPyramid height ctrl xs : - SigSt.toList - (FiltNRG.sumsPosModulatedPyramid - height - (SigSt.fromList SigSt.defaultChunkSize ctrl) - (SigSt.fromList SigSt.defaultChunkSize xs)) : - SigSt.toList - (FiltNRSt.sumsPosModulatedPyramid - height - (SigSt.fromList SigSt.defaultChunkSize ctrl) - (SigSt.fromList SigSt.defaultChunkSize xs)) : - Causal.apply - (FiltNRC.sumsPosModulatedFromPyramid $ - FiltNRSt.pyramid (+) height $ - SigSt.fromList SigSt.defaultChunkSize xs) - ctrl : - [] - -minPosModulated :: - NonNeg.Int -> Sig.T (NonNeg.Int,NonNeg.Int) -> NonEmpty.T Int -> Bool -minPosModulated nheight nctrl xsc = - let ctrl = - map (\(nl,nr) -> - if nl==nr - then (NonNeg.toNumber nl, NonNeg.toNumber nr+1) - else (NonNeg.toNumber nl, NonNeg.toNumber nr)) - nctrl - xs = NonEmpty.toInfiniteList xsc - height = min 10 $ NonNeg.toNumber nheight - in -- trace (show (height, ctrl, xsc)) $ - equalList $ - zipWith FiltNR.minRange (List.tails xs) ctrl : - SigSt.toList - (FiltNRSt.accumulateBinPosModulatedPyramid min height - (SigSt.fromList SigSt.defaultChunkSize ctrl) - (SigSt.fromList SigSt.defaultChunkSize xs)) : - [] - -downSample2 :: - [Int] -> (Int, Sig.T Int) -> Bool -downSample2 lazySize xsc = - let len = Chunky.fromChunks $ map (VP.chunkSize . succ . abs) lazySize - xs = VP.pack len $ cycle $ uncurry (:) xsc - in equalList $ - FiltNRG.downsample2 SigG.defaultLazySize xs : - FiltNRSt.downsample2 xs : - [] - -sumsDownSample2 :: - [Int] -> (Int, Sig.T Int) -> Bool -sumsDownSample2 lazySize xsc = - let len = Chunky.fromChunks $ map (VP.chunkSize . succ . abs) lazySize - xs = VP.pack len $ cycle $ uncurry (:) xsc - in equalList $ - FiltNRG.sumsDownsample2 SigG.defaultLazySize xs : - FiltNRSt.sumsDownsample2 xs : - FiltNRSt.sumsDownsample2Alt xs : - [] - -{- -sumsDownSample2 :: - [VP.ChunkSize] -> (Int, Sig.T Int) -> Bool -sumsDownSample2 lazySize xsc = - let len = Chunky.fromChunks $ filter (0/=) lazySize - xs = VP.pack len $ cycle $ uncurry (:) xsc - in equalList $ - FiltNRG.sumsDownsample2 SigG.defaultLazySize xs : - FiltNRSt.sumsDownsample2 xs : - FiltNRSt.sumsDownsample2Alt xs : - [] --} - -movingAverageModulatedPyramid :: - NonNeg.Int -> Sig.T NonNeg.Int -> - (Stereo.T GF.T, Sig.T (Stereo.T GF.T)) -> Bool -movingAverageModulatedPyramid nheight nctrl xsc = - let ctrl = map NonNeg.toNumber nctrl - xs = uncurry (:) xsc - pack ys = SigSt.fromList SigSt.defaultChunkSize ys - maxC = maximum ctrl - height = min 10 $ NonNeg.toNumber nheight - onegf :: GF.T - onegf = one - in -- trace (show (height, ctrl, xsc)) $ - equalList $ - pack (FiltNR.movingAverageModulatedPyramid onegf - height maxC ctrl (cycle xs)) : - FiltNRG.movingAverageModulatedPyramid onegf - height maxC (pack ctrl) (SigG.cycle $ pack xs) : - FiltNRSt.movingAverageModulatedPyramid onegf - height maxC (pack ctrl) (SigG.cycle $ pack xs) : - [] - - -tests :: [(String, IO ())] -tests = - ("sums", quickCheck sums) : - ("sumRange", quickCheck sumRange) : - ("getRange", quickCheck getRange) : - ("sumsPosModulated", quickCheck sumsPosModulated) : - ("minPosModulated", quickCheck minPosModulated): - ("downSample2", quickCheck downSample2) : - ("sumsDownSample2", quickCheck sumsDownSample2) : - ("movingAverageModulatedPyramid", quickCheck movingAverageModulatedPyramid) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs b/src/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs +++ /dev/null @@ -1,56 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Filter.Allpass (tests) where - -import qualified Synthesizer.Plain.Filter.Recursive.Allpass as Allpass --- import qualified Synthesizer.Plain.Signal as Sig - --- import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty - -import Test.QuickCheck (quickCheck, {- Property, (==>) -}) -import Test.Utility (equalList, ) - --- import qualified Algebra.Module as Module --- import qualified Algebra.RealField as RealField --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive - -import Control.Monad.Trans.State (runState, ) - --- import Debug.Trace (trace, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -{- this will not work due to the poles -parameter :: Double -> Double -> Bool -parameter phase freq = - approxEqual eps phase - (Allpass.makePhase (Allpass.parameter phase freq) freq) --} - - -cascadeStep :: Rational -> Rational -> (Rational, Rational, [Rational]) -> Bool -cascadeStep k u (s0,s1,ns) = - let p = Allpass.Parameter k - s = s0:s1:ns - in equalList $ - runState (Allpass.cascadeStepStack p u) s : - runState (Allpass.cascadeStepRec p u) s : - runState (Allpass.cascadeStepScanl p u) s : - [] - - -cascade :: NonNeg.Int -> Sig.T Rational -> Sig.T Rational -> Bool -cascade order ks xs = - let ps = map Allpass.Parameter ks - n = NonNeg.toNumber order - in Allpass.cascadeState n ps xs == - Allpass.cascadeIterative n ps xs - - -tests :: [(String, IO ())] -tests = - ("cascadeStep", quickCheck cascadeStep) : - ("cascade", quickCheck cascade) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs b/src/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs +++ /dev/null @@ -1,44 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Filter.Hilbert (tests) where - -import qualified Synthesizer.Plain.Filter.Recursive.Hilbert as Hilbert -import qualified Synthesizer.Plain.Filter.Recursive.Allpass as Allpass -import qualified Synthesizer.Plain.Signal as Sig - -import qualified Synthesizer.Causal.Process as Causal - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty - -import Test.QuickCheck (quickCheck, {- Property, (==>) -}) --- import Test.Utility (equalList, ) - --- import qualified Algebra.Module as Module --- import qualified Algebra.RealField as RealField --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive --- import qualified Number.Complex as Complex - -import Data.Tuple.HT (mapPair, ) - --- import Debug.Trace (trace, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -cascade :: NonEmpty.T (Rational, Rational) -> Sig.T Rational -> Bool -cascade ks xs = - let p = uncurry Hilbert.Parameter $ unzip $ - map (mapPair (Allpass.Parameter, Allpass.Parameter)) $ - NonEmpty.toList ks - in Hilbert.run2 p xs == - Causal.apply (Hilbert.causal2 p) xs -{- - in map Complex.real (Hilbert.run2 p xs) == xs --} - - -tests :: [(String, IO ())] -tests = - ("hilbert", quickCheck cascade) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Interpolation.hs b/src/Test/Sound/Synthesizer/Plain/Interpolation.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Interpolation.hs +++ /dev/null @@ -1,343 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Interpolation ( - T, ip, - LinePreserving, lpIp, - tests, - use, useLP, use2, - -- only for debugging - frequencyModulationBackCompare, - frequencyModulationForth0Compare, - frequencyModulationStorableChunkSizeCompare, - frequencyModulationStorableCompare, - ) where - -import qualified Synthesizer.Plain.Interpolation as Interpolation -import qualified Synthesizer.Interpolation.Class as Interpol -import qualified Synthesizer.Interpolation.Custom as ExampleCustom -import qualified Synthesizer.Interpolation.Module as ExampleModule -import qualified Synthesizer.Interpolation as InterpolationCore - -import qualified Synthesizer.Causal.Interpolation as InterpolC -import qualified Synthesizer.Causal.Process as Causal -import qualified Synthesizer.Generic.Filter.NonRecursive as FiltG -import qualified Synthesizer.Generic.Signal as SigG -import qualified Synthesizer.State.Filter.NonRecursive as FiltS -import qualified Synthesizer.State.Signal as SigS - -import qualified Synthesizer.Storable.Filter.NonRecursive as FiltSt -import qualified Synthesizer.Storable.Signal as SigSt - -import Test.QuickCheck (quickCheck, Arbitrary(arbitrary), elements, {- Property, (==>), -} Testable, ) --- import Test.Utility - -import Foreign.Storable (Storable, ) - -import qualified Algebra.VectorSpace as VectorSpace -import qualified Algebra.Module as Module -import qualified Algebra.RealField as RealField -import qualified Algebra.Field as Field -import qualified Algebra.RealRing as RealRing --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty -import qualified Data.List.Match as Match -import Control.Monad (liftM2, ) - -import Test.Utility (equalList, ArbChar, unpackArbString, ) - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - - -instance Arbitrary InterpolationCore.Margin where - arbitrary = - liftM2 InterpolationCore.Margin - (fmap abs arbitrary) - (fmap abs arbitrary) - - -use :: - (Interpolation.T a v -> x) -> - (T a v -> x) -use f ipt = - f (ip ipt) - -useLP :: - (Interpolation.T a v -> x) -> - (LinePreserving a v -> x) -useLP f ipt = - f (lpIp ipt) - -use2 :: - (Interpolation.T a v -> - Interpolation.T a v -> x) -> - (T a v -> - T a v -> x) -use2 f = - use $ \ ipLeap -> - use $ \ ipStep -> - f ipLeap ipStep - - - -data T a v = Cons {name :: String, ip :: Interpolation.T a v} - -instance Show (T a v) where - show x = name x - -instance (Field.C a, Interpol.C a v) => Arbitrary (T a v) where - arbitrary = elements $ - Cons "constant" ExampleCustom.constant : - Cons "linear" ExampleCustom.linear : - Cons "cubic" ExampleCustom.cubic : - [] - - - -data LinePreserving a v = - LPCons {lpName :: String, lpIp :: Interpolation.T a v} - -instance Show (LinePreserving a v) where - show x = lpName x - -instance (Field.C a, Interpol.C a v) => Arbitrary (LinePreserving a v) where - arbitrary = elements $ - LPCons "linear" ExampleCustom.linear : - LPCons "cubic" ExampleCustom.cubic : - [] - - - -constant :: - (Interpol.C a v, Module.C a v, Eq v) => - a -> v -> [v] -> Bool -constant t x0 xs = - equalList $ map ($(x0:xs)) $ map ($t) $ - Interpolation.func ExampleCustom.constant : - Interpolation.func ExampleCustom.piecewiseConstant : - Interpolation.func ExampleModule.constant : - Interpolation.func ExampleModule.piecewiseConstant : - [] - -linear :: - (Interpol.C a v, Module.C a v, Eq v) => - a -> v -> v -> [v] -> Bool -linear t x0 x1 xs = - equalList $ map ($(x0:x1:xs)) $ map ($t) $ - Interpolation.func ExampleCustom.linear : - Interpolation.func ExampleCustom.piecewiseLinear : - Interpolation.func ExampleModule.linear : - Interpolation.func ExampleModule.piecewiseLinear : - [] - -cubic :: - (Interpol.C a v, VectorSpace.C a v, Eq v) => - a -> v -> v -> v -> v -> [v] -> Bool -cubic t x0 x1 x2 x3 xs = - equalList $ map ($(x0:x1:x2:x3:xs)) $ map ($t) $ - Interpolation.func ExampleCustom.cubic : - Interpolation.func ExampleCustom.piecewiseCubic : - Interpolation.func ExampleModule.cubic : - Interpolation.func ExampleModule.cubicAlt : - Interpolation.func ExampleModule.piecewiseCubic : - [] - - -controlAboveOne :: (RealRing.C t) => [t] -> [t] -controlAboveOne = - map ((one+) . abs) - -frequencyModulationForth0 :: - (RealField.C t, Eq v) => - [t] -> [v] -> Bool -frequencyModulationForth0 cs0 xs = - let cs = controlAboveOne cs0 - in Causal.apply - (InterpolC.relative ExampleModule.constant zero - (FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) (SigS.fromList xs))) - (Match.take xs cs) - == Match.take cs xs - -frequencyModulationForth0Compare :: - (RealField.C t, Eq v) => - [t] -> [v] -> ([v], [v], [v]) -frequencyModulationForth0Compare cs0 xs = - let cs = controlAboveOne cs0 - in (Match.take cs - (Causal.apply - (InterpolC.relative ExampleModule.constant zero - (FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) (SigS.fromList xs))) - (Match.take xs cs)), - SigS.toList - (FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) (SigS.fromList xs)), - Match.take cs xs) - - -frequencyModulationForth1 :: - (RealField.C t, Eq v) => - [t] -> [v] -> Bool -frequencyModulationForth1 cs0 xs = - case controlAboveOne cs0 of - [] -> True - (c:cs) -> - Causal.apply - (InterpolC.relative ExampleModule.constant c - (FiltS.inverseFrequencyModulationFloor - (SigS.fromList ((c+one):cs)) (SigS.fromList xs))) - (Match.take xs cs) - == Match.take cs xs - - - -controlBelowOne :: (RealField.C t) => [t] -> [t] -controlBelowOne = - map fraction - - -frequencyModulationBack :: - (RealField.C t, Eq v) => - [t] -> NonEmpty.T v -> Bool -frequencyModulationBack cs0 xs0 = - let cs = controlBelowOne cs0 - xs = NonEmpty.toInfiniteList xs0 - in take (floor (sum cs)) xs == - (SigS.toList $ - FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) - (SigS.fromList $ - Causal.apply - (InterpolC.relative ExampleModule.constant zero - (SigS.fromList xs)) - cs)) - - -frequencyModulationBackCompare :: - (RealField.C t, Eq v) => - [t] -> [v] -> (SigS.T v, SigS.T v) -frequencyModulationBackCompare cs0 xs = - let cs = controlBelowOne cs0 - in (FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) - (SigS.fromList $ - Causal.apply - (InterpolC.relative ExampleModule.constant zero - (SigS.fromList (cycle xs))) - cs), - SigS.fromList $ - Causal.apply - (InterpolC.relative ExampleModule.constant zero - (SigS.fromList (cycle xs))) - cs) - -frequencyModulationGeneric :: - (RealField.C t, Eq v) => - [t] -> [v] -> Bool -frequencyModulationGeneric cs xs = - SigS.toList - (FiltS.inverseFrequencyModulationFloor - (SigS.fromList cs) (SigS.fromList xs)) - == FiltG.inverseFrequencyModulationFloor - SigG.defaultLazySize cs xs - - -makeChunkSize :: Int -> SigSt.ChunkSize -makeChunkSize size = - SigSt.chunkSize (1 + abs size) - -{- -makeExactFraction :: (Int,Int) -> Double -makeExactFraction (n,d) = - fromIntegral n * 2 ^- (- mod (fromIntegral d) 4) --} - -frequencyModulationStorableChunkSize :: - (Storable v, RealField.C t, Eq v) => - Int -> Int -> - Int -> Int -> - [t] -> [v] -> - Bool -frequencyModulationStorableChunkSize size0 size1 xsize0 xsize1 cs xs = - FiltSt.inverseFrequencyModulationFloor - (makeChunkSize size0) cs - (SigSt.fromList (makeChunkSize xsize0) xs) - == - FiltSt.inverseFrequencyModulationFloor - (makeChunkSize size1) cs - (SigSt.fromList (makeChunkSize xsize1) xs) - - -frequencyModulationStorableChunkSizeCompare :: - (Storable v, RealField.C t, Eq v) => - Int -> Int -> - Int -> Int -> - [t] -> [v] -> - (SigSt.T v, SigSt.T v) -frequencyModulationStorableChunkSizeCompare size0 size1 xsize0 xsize1 cs xs = - (FiltSt.inverseFrequencyModulationFloor - (makeChunkSize size0) cs - (SigSt.fromList (makeChunkSize xsize0) xs), - FiltSt.inverseFrequencyModulationFloor - (makeChunkSize size1) cs - (SigSt.fromList (makeChunkSize xsize1) xs)) - - -frequencyModulationStorable :: - (Storable v, RealField.C t, Eq v) => - Int -> Int -> - [t] -> [v] -> - Bool -frequencyModulationStorable size xsize cs xs = - SigSt.toList - (FiltSt.inverseFrequencyModulationFloor (makeChunkSize size) cs - (SigSt.fromList (makeChunkSize xsize) xs)) - == FiltG.inverseFrequencyModulationFloor - SigG.defaultLazySize cs xs - - -frequencyModulationStorableCompare :: - (Storable v, RealField.C t, Eq v) => - Int -> Int -> - [t] -> [v] -> - ([v], SigSt.T v) -frequencyModulationStorableCompare size xsize cs xs = - (FiltG.inverseFrequencyModulationFloor - SigG.defaultLazySize cs xs, - FiltSt.inverseFrequencyModulationFloor (makeChunkSize size) cs - (SigSt.fromList (makeChunkSize xsize) xs)) - - - -testRational :: - (Testable t) => - (Rational -> Rational -> t) -> IO () -testRational = quickCheck - -testFM :: - (Testable t, Arbitrary (sigX ArbChar), Show (sigX ArbChar)) => - ([Rational] -> sigX ArbChar -> t) -> IO () -testFM = quickCheck - -tests :: [(String, IO ())] -tests = - ("constant", testRational constant) : - ("linear", testRational linear ) : - ("cubic", testRational cubic ) : - ("frequencyModulationForth0", testFM frequencyModulationForth0) : - ("frequencyModulationForth1", testFM frequencyModulationForth1) : - ("frequencyModulationBack", testFM frequencyModulationBack) : - ("frequencyModulationGeneric", testFM frequencyModulationGeneric) : - ("frequencyModulationStorableChunkSize", - quickCheck (\size0 size1 xsize0 xsize1 cs xs -> - frequencyModulationStorableChunkSize size0 size1 xsize0 xsize1 - (cs::[Rational]) (unpackArbString xs))) : - ("frequencyModulationStorable", - quickCheck (\size xsize cs xs -> - frequencyModulationStorable size xsize - (cs::[Rational]) (unpackArbString xs))) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/NonEmpty.hs b/src/Test/Sound/Synthesizer/Plain/NonEmpty.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/NonEmpty.hs +++ /dev/null @@ -1,34 +0,0 @@ -module Test.Sound.Synthesizer.Plain.NonEmpty where - -import Test.QuickCheck (Arbitrary, arbitrary, ) -import Control.Monad (liftM2, ) - - -data T a = Cons a [a] - -toList :: T a -> [a] -toList (Cons x xs) = - (x:xs) - -toInfiniteList :: T a -> [a] -toInfiniteList = - cycle . toList - -instance Functor T where - fmap f (Cons x xs) = - Cons (f x) (map f xs) - -instance Arbitrary a => Arbitrary (T a) where - arbitrary = liftM2 Cons arbitrary arbitrary - -instance Show a => Show (T a) where - showsPrec p (Cons x xs) = - showsPrec p (x:xs) - -{- -instance Show a => Show (T a) where - showsPrec p (Cons x xs) = - showParen (p >= 10) $ - showString "cycle " . - showsPrec 11 (x:xs) --} diff --git a/src/Test/Sound/Synthesizer/Plain/Oscillator.hs b/src/Test/Sound/Synthesizer/Plain/Oscillator.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Oscillator.hs +++ /dev/null @@ -1,39 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Oscillator (tests) where - -import qualified Synthesizer.Plain.Oscillator as Osci -import qualified Synthesizer.Basic.Wave as Wave --- import qualified Synthesizer.Plain.Interpolation as Interpolation - -import qualified Test.Sound.Synthesizer.Plain.Wave as WaveTest --- import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest - -import Test.QuickCheck (quickCheck, {- Property, (==>), -} ) - -import qualified Algebra.RealField as RealField - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - - -phaseShapeMod :: (RealField.C a, Eq b) => (Wave.T a b) -> a -> [a] -> Bool -phaseShapeMod wave freq phases = - Osci.phaseMod wave freq phases == - Osci.shapeMod (Wave.phaseOffset wave) zero freq phases - -phaseShapeModRational :: - WaveTest.Ring Rational -> Integer -> Integer -> [Integer] -> Bool -phaseShapeModRational w denom0 freq0 phases0 = - let denom = 1 + abs denom0 - freq = freq0 % denom - phases = map (% denom) phases0 - in phaseShapeMod (WaveTest.ringWave w) freq phases - - - -tests :: [(String, IO ())] -tests = - ("phaseShapeModRational", quickCheck phaseShapeModRational) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/ToneModulation.hs b/src/Test/Sound/Synthesizer/Plain/ToneModulation.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/ToneModulation.hs +++ /dev/null @@ -1,478 +0,0 @@ -module Test.Sound.Synthesizer.Plain.ToneModulation (tests, ) where - -import Test.Sound.Synthesizer.Basic.ToneModulation ( - minLength, - minLengthMargin, - shapeLimits, - testRationalLineIp, - testRationalIp, - ) - -import qualified Synthesizer.Plain.Oscillator as Osci -import qualified Synthesizer.Plain.Interpolation as Interpolation -import qualified Synthesizer.Plain.ToneModulation as ToneModL -import qualified Synthesizer.Plain.Wave as WaveL -import Synthesizer.Interpolation (marginNumber, ) - -import qualified Synthesizer.Basic.Wave as Wave -import qualified Synthesizer.Basic.Phase as Phase - -import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty -import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest - -import Test.QuickCheck (quickCheck, Property, (==>), ) -import Test.Utility (ArbChar, ) - -import qualified Number.NonNegative as NonNeg -import qualified Number.NonNegativeChunky as Chunky - -import qualified Algebra.RealTranscendental as RealTrans -import qualified Algebra.Module as Module -import qualified Algebra.RealField as RealField -import qualified Algebra.Additive as Additive -import qualified Algebra.ZeroTestable as ZeroTestable - -import Data.List.HT (isAscending, ) -import Data.Ord.HT (limit, ) -import Data.Tuple.HT (mapPair, mapSnd, ) -import qualified Data.List as List - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -{- -Properties that do not hold: - commutativity of limitRelativeShapes and integrateFractional: - Does not hold because when you clip the integral skips at the end, - you would have to clear the fractional part, too. --} - - - -absolutize :: (Additive.C a) => a -> [a] -> [a] -absolutize = scanl (+) - -limitMinRelativeValues :: - Int -> Int -> [NonNeg.Int] -> Bool -limitMinRelativeValues xMin x0 xsnn = - let xs = map NonNeg.toNumber xsnn - in map (max xMin) (absolutize x0 xs) == - uncurry absolutize (ToneModL.limitMinRelativeValues xMin x0 xs) - -limitMaxRelativeValues :: - Int -> Int -> [NonNeg.Int] -> Bool -limitMaxRelativeValues xMax x0 xsnn = - let xs = map NonNeg.toNumber xsnn - in map (min xMax) (absolutize x0 xs) == - uncurry absolutize (ToneModL.limitMaxRelativeValues xMax x0 xs) - -limitMaxRelativeValuesNonNeg :: - Int -> Int -> [NonNeg.Int] -> Bool -limitMaxRelativeValuesNonNeg xMax x0 xsnn = - let xs = map NonNeg.toNumber xsnn - in map (min xMax) (absolutize x0 xs) == - uncurry absolutize (ToneModL.limitMaxRelativeValuesNonNeg xMax x0 xs) - --- chunky type is not necessary here but testing it a little is not wrong -limitMinRelativeValuesIdentity :: - Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool -limitMinRelativeValuesIdentity x0 xs = - (x0,xs) == ToneModL.limitMinRelativeValues 0 x0 xs - -limitMaxRelativeValuesIdentity :: - Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool -limitMaxRelativeValuesIdentity x0 xs = - let inf = 1 + inf - in (x0,xs) == ToneModL.limitMaxRelativeValues inf x0 xs - -limitMaxRelativeValuesNonNegIdentity :: - Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool -limitMaxRelativeValuesNonNegIdentity x0 xs = - let inf = 1 + inf - in (x0,xs) == ToneModL.limitMaxRelativeValuesNonNeg inf x0 xs - -limitMaxRelativeValuesInfinity :: - Chunky.T NonNeg.Int -> NonEmpty.T (Chunky.T NonNeg.Int) -> Bool -limitMaxRelativeValuesInfinity x0 ixs = - let inf = 1 + inf - ys = NonEmpty.toInfiniteList ixs - (z0,zs) = ToneModL.limitMaxRelativeValues inf x0 ys - in (x0, take 100 ys) == (z0, take 100 zs) - -limitMaxRelativeValuesNonNegInfinity :: - Chunky.T NonNeg.Int -> NonEmpty.T (Chunky.T NonNeg.Int) -> Bool -limitMaxRelativeValuesNonNegInfinity x0 ixs = - let inf = 1 + inf - ys = NonEmpty.toInfiniteList ixs - (z0,zs) = ToneModL.limitMaxRelativeValuesNonNeg inf x0 ys - in (x0, take 100 ys) == (z0, take 100 zs) - - -dropRem :: Eq a => NonNeg.Int -> [a] -> Bool -dropRem nn xs = - let n = NonNeg.toNumber nn - in map (flip ToneModL.dropRem xs) [0 .. n + length xs] == - map ((,) 0) (List.tails xs) ++ map (flip (,) []) [1..n] - - -sampledToneSine :: (RealTrans.C a, Module.C a a) => - NonNeg.T a -> NonNeg.Int -> a -> a -> a -> Bool -sampledToneSine periodNN ext phase0 shape phase = - let ipLeap = Interpolation.cubic - ipStep = Interpolation.cubic - ten = fromInteger 10 - period = ten + NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (Osci.staticSine phase0 (recip period)) - in abs (WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` (Phase.fromRepresentative phase) - - head (Osci.staticSine (phase0+phase) zero)) < ten ^- (-2) - - -sampledToneSineList :: (RealTrans.C a, Module.C a a) => - NonNeg.T a -> NonNeg.Int -> a -> a -> [a] -> [a] -> Bool -sampledToneSineList periodNN ext origPhase phase shapes freqs = - let ipLeap = Interpolation.cubic - ipStep = Interpolation.cubic - ten = fromInteger 10 - period = ten + NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (Osci.staticSine origPhase (recip period)) - in all ((< ten ^- (-2)) . abs) $ - zipWith (-) - (Osci.shapeFreqMod (WaveL.sampledTone ipLeap ipStep period tone) - phase shapes freqs) - (Osci.freqModSine (origPhase+phase) freqs) - - -sampledToneLinear :: (RealField.C a, Module.C a v, Eq v) => - InterpolationTest.LinePreserving a v -> - InterpolationTest.LinePreserving a v -> - NonNeg.T a -> NonNeg.Int -> (v,v) -> a -> Phase.T a -> Property -sampledToneLinear = - InterpolationTest.useLP $ \ ipLeap -> - InterpolationTest.useLP $ \ ipStep -> - \ periodNN ext (i,d) shape phase -> - let period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - ramp = take len (List.iterate (d+) i) - limits = - mapPair (fromIntegral, fromIntegral) $ - shapeLimits ipLeap ipStep periodInt len - in period /= zero ==> - -- should be (fraction phase), right? - WaveL.sampledTone ipLeap ipStep period ramp shape `Wave.apply` phase == - i + limit limits shape *> d -{- -let len=100; period=1/0.06::Double; ip = Interpolation.linear in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (0,fromIntegral len)) [\s -> WaveL.sampledTone ip ip period (take len $ iterate (1+) (0::Double)) s 0, limit (mapPair (fromIntegral, fromIntegral) $ shapeLimits ip ip (round period::Int) len)] --} - -sampledToneStair :: (RealField.C a, Module.C a v, Eq v) => - InterpolationTest.LinePreserving a v -> - NonNeg.Int -> NonNeg.Int -> (v,v) -> a -> Property -sampledToneStair = - InterpolationTest.useLP $ \ ipLeap - periodIntNN ext (i,d) shape -> - let ipStep = Interpolation.constant - periodInt = NonNeg.toNumber periodIntNN - period = fromIntegral periodInt - len0 = minLength ipLeap ipStep periodInt ext - (rep,rm) = divMod (negate len0) periodInt - len = len0 + rm - stair = - concatMap (replicate periodInt) $ - take (negate rep) (List.iterate (period*>d+) i) - limits = - mapPair (fromIntegral, fromIntegral) $ - shapeLimits ipLeap ipStep periodInt len - in periodInt /= zero ==> - WaveL.sampledTone ipLeap ipStep period stair shape `Wave.apply` zero == - i + limit limits shape *> d -{- -let len=periodInt*rep; rep=10; periodInt = 14::Int; period=fromIntegral periodInt; ipl = Interpolation.linear; ipc = Interpolation.constant in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (-10,10+fromIntegral len)) [\s -> WaveL.sampledTone ipl ipc period (concatMap (replicate periodInt) $ take rep $ iterate (period+) (0::Double)) s 0, limit (mapPair (fromIntegral, fromIntegral) $ shapeLimits ipl ipc periodInt len)] --} - -{- -sampledToneSaw :: (RealField.C a, Module.C a v, Eq v) => - InterpolationTest.LinePreserving a v -> - InterpolationTest.T a v -> - NonNeg.Int -> NonNeg.Int -> (v,v) -> a -> a -> Property -sampledToneSaw iptLeap iptStep periodIntNN ext (i,d) shape phase = - let ipLeap = InterpolationTest.lpIp iptLeap - ipStep = InterpolationTest.ip iptStep - periodInt = NonNeg.toNumber periodIntNN - period = fromIntegral periodInt - len0 = minLength ipLeap ipStep periodInt ext - rep = negate $ div (negate len0) periodInt - saw = - concat $ replicate rep $ - take periodInt $ List.iterate (d+) i - in periodInt /= zero ==> - WaveL.sampledTone ipLeap ipStep period saw shape phase == - i + fraction phase *> d --} - -sampledToneStatic :: (RealField.C a, Eq v) => - InterpolationTest.T a v -> - InterpolationTest.T a v -> - NonNeg.Int -> (v,[v]) -> a -> a -> Property -sampledToneStatic = - InterpolationTest.use2 $ \ ipLeap ipStep - ext (x,xs) shape phase -> - let wave = x:xs - periodInt = length wave - period = fromIntegral periodInt - len = minLength ipLeap ipStep periodInt ext - rep = negate $ div (negate len) periodInt - tone = concat $ replicate rep wave - in period /= zero ==> - WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` (Phase.fromRepresentative phase) == - Interpolation.cyclicPad Interpolation.single ipStep (phase*period) wave -{- -let wave = [1,-1,0.5,-0.5::Double]; period = fromIntegral (length wave) :: Double; ip = Interpolation.linear in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (-1,3)) [WaveL.sampledTone ip ip period (concat $ replicate 3 wave) 0.3, \phase -> Interpolation.cyclicPad Interpolation.single Interpolation.linear (phase*period) wave] --} - - - -shapeFreqModFromSampledToneLimitIdentity :: (RealField.C t) => - Interpolation.Margin -> - Interpolation.Margin -> - NonNeg.Int -> NonEmpty.T y -> (t, NonEmpty.T (NonNeg.T t)) -> Bool -shapeFreqModFromSampledToneLimitIdentity - marginLeap marginStep periodIntNN ixs (shape0,shapesNN) = - let periodInt = NonNeg.toNumber periodIntNN - shapes = fmap NonNeg.toNumber shapesNN - a = snd - (ToneModL.limitRelativeShapes - marginLeap marginStep - periodInt (NonEmpty.toInfiniteList ixs) - (shape0, NonEmpty.toInfiniteList shapes)) !! 100 - in a == a - - -oscillatorCoords :: (RealField.C t) => - NonNeg.Int -> NonNeg.T t -> t -> Phase.T t -> [NonNeg.T t] -> [t] -> Property -oscillatorCoords - periodIntNN periodNN shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - periodRound = fromIntegral periodInt - coords = - ToneModL.oscillatorCoords - periodInt period - (shape0, shapes) (phase, freqs) - in period /= zero && periodInt /= zero ==> - all - (\(skip,(k,(qShape,qWave))) -> - skip >= zero && - isAscending [negate periodInt, k, zero] && - isAscending [zero, qShape, one] && - isAscending [zero, qWave, periodRound]) - (tail coords) - - -shapeFreqModFromSampledToneCoordsIdentity :: - (RealField.C t, ZeroTestable.C t) => - NonNeg.Int -> NonNeg.T t -> (t, [NonNeg.T t]) -> Property -shapeFreqModFromSampledToneCoordsIdentity - periodIntNN periodNN (shape0,shapesNN) = - let period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - shapes = map NonNeg.toNumber shapesNN - phase = Phase.fromRepresentative $ shape0 / period - freqs = map (/period) shapes - in period /= zero ==> - all - (isZero . fst . snd . snd) - (ToneModL.oscillatorCoords - periodInt period (shape0, shapes) (phase, freqs)) - - -shapeFreqModFromSampledTone :: (RealField.C t, Eq v) => - InterpolationTest.T t v -> - InterpolationTest.T t v -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> - Property -shapeFreqModFromSampledTone = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs shape0 phase shapesNN freqs -> - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledToneA = - Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone - shape0 phase shapes freqs - resampledToneB = - Osci.shapeFreqMod - (WaveL.sampledTone ipLeap ipStep period tone) - phase (scanl (+) shape0 shapes) freqs - in period /= zero ==> - resampledToneA == resampledToneB -{- -let len=100; period=1/0.06::Double; ip = Interpolation.linear; tone = take len $ iterate (1+) (0::Double); shape0=0; shapes = replicate 100 1; in GNUPlot.plotLists [] [Osci.shapeFreqMod (WaveL.sampledTone ip ip period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ip ip period tone shape0 0 shapes (repeat 0)] -*Test.Sound.Synthesizer.Plain.Oscillator> let len=100; period=1/0.06::Double; ip = Interpolation.linear; tone = take len $ iterate (1+) (0::Double); shape0=0; shapes = concat $ replicate 50 [1.5,0.5]; in GNUPlot.plotLists [] [Osci.shapeFreqMod (WaveL.sampledTone ip ip period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ip ip period tone shape0 0 shapes (repeat 0)] -*Test.Sound.Synthesizer.Plain.Oscillator> let len=100; period=1/0.06::Rational; ipLeap = Interpolation.linear; ipStep = Interpolation.constant; tone = take len $ iterate (1+) (0::Rational); shape0=0; shapes = concat $ replicate 50 [1.5,0.5]; in GNUPlot.plotLists [] (map (map (\x -> fromRational' x :: Double)) [Osci.shapeFreqMod (WaveL.sampledTone ipLeap ipStep period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone shape0 0 shapes (repeat 0)]) --} - - -shapePhaseFreqModFromSampledTone :: (RealField.C t, Eq v) => - InterpolationTest.T t v -> - InterpolationTest.T t v -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> [t] -> - Property -shapePhaseFreqModFromSampledTone = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs shape0 phase shapesNN phaseDistorts freqs -> - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - resampledToneA = - Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone - shape0 phase shapes phaseDistorts freqs - resampledToneB = - Osci.shapeFreqMod - (uncurry $ - Wave.phaseOffset . - WaveL.sampledTone ipLeap ipStep period tone) - phase (zip (scanl (+) shape0 shapes) phaseDistorts) freqs - in period /= zero ==> - resampledToneA == resampledToneB - - -oscillatorCells :: (RealField.C t, Eq v) => - Interpolation.Margin -> - Interpolation.Margin -> - NonNeg.Int -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - t -> t -> [NonNeg.T t] -> [t] -> - Property -oscillatorCells - marginLeap marginStep periodIntNN periodNN ext ixs shape0 phase shapesNN freqs = - let shapes = map NonNeg.toNumber shapesNN - period = NonNeg.toNumber periodNN - periodInt = NonNeg.toNumber periodIntNN - len = minLengthMargin marginLeap marginStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - crop = cropCell marginLeap marginStep - resampledToneA = - ToneModL.oscillatorCells - marginLeap marginStep periodInt period tone - (shape0, shapes) (Phase.fromRepresentative phase, freqs) - resampledToneB = - Osci.shapeFreqMod - (Wave.Cons . ToneModL.sampledToneCell - (ToneModL.makePrototype marginLeap marginStep - periodInt period tone)) - phase (scanl (+) shape0 shapes) freqs - in period /= zero && - periodInt /= zero && - marginNumber marginLeap > zero && - marginNumber marginStep > zero ==> - map crop resampledToneA == map crop resampledToneB - -cropCell :: - Interpolation.Margin -> - Interpolation.Margin -> - ((t,t), ToneModL.Cell v) -> ((t,t), ToneModL.Cell v) -cropCell ipLeap ipStep = - mapSnd - (take (marginNumber ipStep) . - map (take (marginNumber ipLeap))) - - -shapeFreqModFromSampledToneIdentity :: (RealField.C t, Eq v) => - InterpolationTest.T t v -> - InterpolationTest.T t v -> - NonNeg.T t -> - NonNeg.Int -> NonEmpty.T v -> - Property -shapeFreqModFromSampledToneIdentity = - InterpolationTest.use2 $ \ ipLeap ipStep - periodNN ext ixs -> - let period = NonNeg.toNumber periodNN - periodInt = round period - len = minLength ipLeap ipStep periodInt ext - tone = take len (NonEmpty.toInfiniteList ixs) - shape0 = zero - shapes = repeat one - phase = zero - freqs = repeat (recip period) - (n0,n1) = - shapeLimits ipLeap ipStep periodInt len - - resampledTone = - Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone - shape0 phase shapes freqs - in period /= zero ==> - and (drop n0 (take (succ n1) (zipWith (==) resampledTone tone))) - - -tests :: [(String, IO ())] -tests = - ("limitMinRelativeValues", quickCheck limitMinRelativeValues) : - ("limitMaxRelativeValues", quickCheck limitMaxRelativeValues) : - ("limitMaxRelativeValuesNonNeg", - quickCheck limitMaxRelativeValuesNonNeg) : - ("limitMinRelativeValuesIdentity", - quickCheck limitMinRelativeValuesIdentity) : - ("limitMaxRelativeValuesIdentity", - quickCheck limitMaxRelativeValuesIdentity) : - ("limitMaxRelativeValuesNonNegIdentity", - quickCheck limitMaxRelativeValuesNonNegIdentity) : - ("limitMaxRelativeValuesInfinity", - quickCheck limitMaxRelativeValuesInfinity) : - ("limitMaxRelativeValuesNonNegInfinity", - quickCheck limitMaxRelativeValuesNonNegInfinity) : - ("dropRem", quickCheck (dropRem :: NonNeg.Int -> [ArbChar] -> Bool)) : - ("sampledToneSine", - quickCheck (\period -> sampledToneSine (period :: NonNeg.Double))) : - ("sampledToneSineList", - quickCheck (\period -> sampledToneSineList (period :: NonNeg.Double))) : - ("sampledToneLinear", - testRationalLineIp sampledToneLinear) : - ("sampledToneStair", - testRationalLineIp sampledToneStair) : -{- - ("sampledToneSaw", - testRationalLineIp sampledToneSaw) : --} - ("sampledToneStatic", - testRationalIp sampledToneStatic) : - ("shapeFreqModFromSampledToneLimitIdentity", - quickCheck (\ml ms p ixs (t,ts) -> - shapeFreqModFromSampledToneLimitIdentity ml ms p - (ixs::NonEmpty.T Rational) (t::Rational,ts))) : - ("oscillatorCoords", - quickCheck (\periodInt period -> - oscillatorCoords - periodInt (period :: NonNeg.Rational))) : - ("shapeFreqModFromSampledToneCoordsIdentity", - quickCheck (\periodInt period -> - shapeFreqModFromSampledToneCoordsIdentity - periodInt (period :: NonNeg.Rational))) : - ("shapeFreqModFromSampledTone", - testRationalIp shapeFreqModFromSampledTone) : - ("shapePhaseFreqModFromSampledTone", - testRationalIp shapePhaseFreqModFromSampledTone) : - ("oscillatorCells", - quickCheck (\ml ms periodInt period ext ixs -> - oscillatorCells ml ms periodInt (period :: NonNeg.Rational) - ext (ixs :: NonEmpty.T ArbChar))) : - ("shapeFreqModFromSampledToneIdentity", - testRationalIp shapeFreqModFromSampledToneIdentity) : - [] diff --git a/src/Test/Sound/Synthesizer/Plain/Wave.hs b/src/Test/Sound/Synthesizer/Plain/Wave.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Plain/Wave.hs +++ /dev/null @@ -1,75 +0,0 @@ -module Test.Sound.Synthesizer.Plain.Wave (Ring, ringWave, tests) where - -import qualified Synthesizer.Basic.Wave as Wave -import qualified Synthesizer.Basic.Phase as Phase - -import Test.QuickCheck (quickCheck, Arbitrary(arbitrary), elements, oneof, choose, {- Property, (==>), -} ) --- import Test.Utility - -import qualified Number.NonNegative as NonNeg - -import qualified Algebra.RealTranscendental as RealTrans -import qualified Algebra.Ring as Ring - -import Control.Monad (liftM, liftM2, ) -import System.Random (Random) - - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - - - -data Ring a = Ring {ringName :: String, ringWave :: Wave.T a a} - -instance Show (Ring a) where - show = ringName - -instance (Ord a, Ring.C a) => Arbitrary (Ring a) where - arbitrary = elements $ - Ring "saw" Wave.saw : - Ring "square" Wave.square : - Ring "triangle" Wave.triangle : - [] - - - - -data ZeroDCOffset a = ZeroDCOffset {zdcName :: String, zdcWave :: Wave.T a a} - -instance Show (ZeroDCOffset a) where - show = zdcName - -instance (RealTrans.C a, Random a) => Arbitrary (ZeroDCOffset a) where - arbitrary = - let cons n w = return (ZeroDCOffset n w) - in oneof $ - cons "sine" Wave.sine : - cons "saw" Wave.saw : - cons "square" Wave.square : - cons "triangle" Wave.triangle : - liftM - (ZeroDCOffset "squareBalanced" . Wave.squareBalanced) - (choose (negate one, one)) : - liftM2 - (\w r -> ZeroDCOffset "trapezoidBalanced" (Wave.trapezoidBalanced w r)) - (choose (zero, one)) - (choose (negate one, one)) : - [] - - -zeroDCOffset :: ZeroDCOffset Double -> NonNeg.Int -> Bool -zeroDCOffset w periodIntNN = - let periodInt = 100 + NonNeg.toNumber periodIntNN - period = fromIntegral periodInt - xs = take periodInt $ map Phase.fromRepresentative $ - map (/period) $ iterate (1+) 0.5 - in abs (sum (map (Wave.apply (zdcWave w)) xs)) < period / fromInteger 100 - - -tests :: [(String, IO ())] -tests = - ("zeroDCOffset", quickCheck zeroDCOffset) : - [] diff --git a/src/Test/Sound/Synthesizer/Storable/Cut.hs b/src/Test/Sound/Synthesizer/Storable/Cut.hs deleted file mode 100644 --- a/src/Test/Sound/Synthesizer/Storable/Cut.hs +++ /dev/null @@ -1,40 +0,0 @@ -module Test.Sound.Synthesizer.Storable.Cut (tests) where - -import qualified Synthesizer.Storable.Cut as CutSt -import qualified Synthesizer.Storable.Signal as SigSt - -import qualified Synthesizer.Plain.Cut as Cut -import qualified Synthesizer.Plain.Signal as Sig - -import qualified Data.EventList.Relative.TimeBody as EventList - --- import qualified Algebra.RealRing as RealRing --- import qualified Algebra.Ring as Ring --- import qualified Algebra.Additive as Additive - -import qualified Number.NonNegative as NonNeg - -import Test.QuickCheck (quickCheck, ) -import Test.Utility (equalList, ) - -import NumericPrelude.Numeric -import NumericPrelude.Base -import Prelude () - - -arrange :: NonNeg.Int -> EventList.T NonNeg.Int (Sig.T Int) -> Bool -arrange nnChunkSize evs = - let chunkSize = SigSt.chunkSize $ 1 + NonNeg.toNumber nnChunkSize - sevs = EventList.mapBody (SigSt.fromList chunkSize) evs - in equalList $ - SigSt.fromList chunkSize (Cut.arrange evs) : - CutSt.arrangeAdaptive chunkSize sevs : - CutSt.arrangeList chunkSize sevs : - CutSt.arrangeEquidist chunkSize sevs : - [] - - -tests :: [(String, IO ())] -tests = - ("arrange", quickCheck arrange) : - [] diff --git a/src/Test/Utility.hs b/src/Test/Utility.hs deleted file mode 100644 --- a/src/Test/Utility.hs +++ /dev/null @@ -1,69 +0,0 @@ -{-# LANGUAGE NoImplicitPrelude #-} -module Test.Utility where - -import Test.QuickCheck (Arbitrary(arbitrary)) - -import qualified Number.Complex as Complex - -import qualified Algebra.RealRing as RealRing - -import qualified Data.List.HT as ListHT -import qualified Data.Char as Char - -import NumericPrelude.Base -import NumericPrelude.Numeric - - -equalList :: Eq a => [a] -> Bool -equalList xs = - and (ListHT.mapAdjacent (==) xs) - - -approxEqual :: (RealRing.C a) => a -> a -> a -> Bool -approxEqual eps x y = - 2 * abs (x-y) <= eps * (abs x + abs y) - -approxEqualAbs :: (RealRing.C a) => a -> a -> a -> Bool -approxEqualAbs eps x y = - abs (x-y) <= eps - -approxEqualListRel :: (RealRing.C a) => a -> [a] -> Bool -approxEqualListRel eps xs = - let n = fromIntegral $ length xs - in approxEqualListAbs (eps * n * sum (map abs xs)) xs - -approxEqualListAbs :: (RealRing.C a) => a -> [a] -> Bool -approxEqualListAbs eps xs = - let n = fromIntegral $ length xs - s = sum xs - in sum (map (\x -> abs (n*x-s)) xs) <= eps - - -approxEqualComplex :: - (RealRing.C a) => - a -> Complex.T a -> Complex.T a -> Bool -approxEqualComplex eps x y = - 2 * Complex.magnitudeSqr (x-y) - <= eps^2 * (Complex.magnitudeSqr x + Complex.magnitudeSqr y) - -approxEqualComplexAbs :: - (RealRing.C a) => - a -> Complex.T a -> Complex.T a -> Bool -approxEqualComplexAbs eps x y = - Complex.magnitudeSqr (x-y) <= eps^2 - - --- see event-list - -newtype ArbChar = ArbChar Char - deriving (Eq, Ord) - -instance Show ArbChar where - showsPrec n (ArbChar c) = showsPrec n c - -instance Arbitrary ArbChar where - arbitrary = fmap (ArbChar . Char.chr . (32+) . flip mod 96) arbitrary - -unpackArbString :: [ArbChar] -> String -unpackArbString = - map (\(ArbChar c) -> c) diff --git a/synthesizer-core.cabal b/synthesizer-core.cabal --- a/synthesizer-core.cabal +++ b/synthesizer-core.cabal @@ -1,5 +1,5 @@ Name: synthesizer-core -Version: 0.7 +Version: 0.7.0.1 License: GPL License-File: LICENSE Author: Henning Thielemann @@ -25,7 +25,7 @@ Stability: Experimental Tested-With: GHC==6.4.1, GHC==6.8.2, GHC==6.10.4, GHC==6.12.3 Tested-With: GHC==7.0.4, GHC==7.2.1, GHC==7.4.2, GHC==7.6.3 -Cabal-Version: >=1.6 +Cabal-Version: >=1.14 Build-Type: Simple Extra-Source-Files: @@ -48,7 +48,7 @@ Source-Repository this - Tag: 0.7 + Tag: 0.7.0.1 Type: darcs Location: http://code.haskell.org/synthesizer/core/ @@ -87,10 +87,11 @@ -- also warns about NumericPrelude import: -fwarn-missing-import-lists GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP GHC-Options: -Wall - Hs-source-dirs: src + Hs-source-dirs: src, private Exposed-modules: Synthesizer.Storage @@ -230,15 +231,30 @@ Executable test - If !flag(buildTests) + If flag(buildTests) + Build-Depends: + synthesizer-core, + storablevector, + storable-tuple, + event-list, + non-empty, + non-negative, + utility-ht, + numeric-prelude, + QuickCheck, + random, + containers, + base + Else Buildable: False GHC-Options: -Wall -fwarn-tabs -fwarn-incomplete-record-updates - Hs-Source-Dirs: src + Hs-Source-Dirs: test, private If impl(ghc>=7.0) GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP Other-Modules: Test.Utility @@ -268,10 +284,12 @@ Executable fouriertest If flag(buildProfilers) Build-Depends: + synthesizer-core, + numeric-prelude, + timeit >=1.0 && <1.1, storablevector >=0.2.7 && <0.3, - utility-ht >=0.0.5 && <0.1, storable-tuple >=0.0.1 && <0.1, - timeit >=1.0 && <1.1, + utility-ht >=0.0.5 && <0.1, base >=4 && <5 Else Buildable: False @@ -279,56 +297,85 @@ If impl(ghc>=7.0) GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP GHC-Options: -Wall GHC-Prof-Options: -auto-all - Hs-Source-Dirs: speedtest, src + Hs-Source-Dirs: speedtest Main-Is: Fourier.hs Executable speedtest - If !flag(buildProfilers) + If flag(buildProfilers) + Build-Depends: + synthesizer-core, + numeric-prelude, + old-time, + directory, + binary, + bytestring, + utility-ht, + base + Else Buildable: False If impl(ghc>=7.0) GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP GHC-Options: -Wall -fexcess-precision If flag(optimizeAdvanced) GHC-Options: -optc-ffast-math -optc-O3 -- -funfolding-use-threshold=20 -funfolding-creation-threshold=100 -- -optc-march=pentium4 -optc-mfpmath=sse - Hs-Source-Dirs: speedtest, src + Hs-Source-Dirs: speedtest Main-Is: SpeedTest.hs Executable speedtest-exp - If !flag(buildProfilers) + If flag(buildProfilers) + Build-Depends: + synthesizer-core, + storablevector, + binary, + bytestring, + array, + base + If flag(splitBase) + Build-Depends: + old-time >= 1.0 && < 1.2, + directory >= 1.0 && < 1.3 + Else Buildable: False If impl(ghc>=7.0) GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP GHC-Options: -Wall -fexcess-precision - Hs-Source-Dirs: speedtest, src + Hs-Source-Dirs: speedtest Main-Is: SpeedTestExp.hs - If flag(splitBase) - Build-Depends: - old-time >= 1.0 && < 1.2, - directory >= 1.0 && < 1.3 Executable speedtest-simple - If !flag(buildProfilers) + If flag(buildProfilers) + Build-Depends: + synthesizer-core, + binary, + bytestring, + old-time, + base + Else Buildable: False If impl(ghc>=7.0) GHC-Options: -fwarn-unused-do-bind CPP-Options: -DNoImplicitPrelude=RebindableSyntax - Extensions: CPP + Default-Language: Haskell2010 + Default-Extensions: CPP GHC-Options: -Wall - Hs-Source-Dirs: speedtest, src + Hs-Source-Dirs: speedtest Main-Is: SpeedTestSimple.hs diff --git a/test/Test/Main.hs b/test/Test/Main.hs new file mode 100644 --- /dev/null +++ b/test/Test/Main.hs @@ -0,0 +1,49 @@ +module Main where + +import qualified Test.Sound.Synthesizer.Plain.Analysis as Analysis +import qualified Test.Sound.Synthesizer.Plain.Control as Control +import qualified Test.Sound.Synthesizer.Plain.Filter as Filter +import qualified Test.Sound.Synthesizer.Plain.Interpolation as Interpolation +import qualified Test.Sound.Synthesizer.Plain.Oscillator as Oscillator +import qualified Test.Sound.Synthesizer.Plain.Wave as Wave +import qualified Test.Sound.Synthesizer.Basic.NumberTheory as NumberTheory +import qualified Test.Sound.Synthesizer.Basic.ToneModulation as ToneModulation +import qualified Test.Sound.Synthesizer.Plain.ToneModulation as ToneModulationL +import qualified Test.Sound.Synthesizer.Generic.ToneModulation as ToneModulationG +import qualified Test.Sound.Synthesizer.Generic.Permutation as Permutation +import qualified Test.Sound.Synthesizer.Generic.Fourier as Fourier +import qualified Test.Sound.Synthesizer.Generic.FourierInteger as FourierInteger +import qualified Test.Sound.Synthesizer.Generic.Filter as FilterG +import qualified Test.Sound.Synthesizer.Generic.Cut as CutG +import qualified Test.Sound.Synthesizer.Causal.Analysis as AnalysisC +import qualified Test.Sound.Synthesizer.Storable.Cut as Cut + +import Data.Tuple.HT (mapFst, ) + + +prefix :: String -> [(String, IO ())] -> [(String, IO ())] +prefix msg = + map (mapFst (\str -> msg ++ "." ++ str)) + +main :: IO () +main = + mapM_ (\(msg,io) -> putStr (msg++": ") >> io) $ + concat $ + prefix "Plain.Analysis" Analysis.tests : + prefix "Plain.Control" Control.tests : + prefix "Plain.Filter" Filter.tests : + prefix "Plain.Interpolation" Interpolation.tests : + prefix "Plain.Oscillator" Oscillator.tests : + prefix "Plain.Wave" Wave.tests : + prefix "Storable.Cut" Cut.tests : + prefix "Generic.Cut" CutG.tests : + prefix "Basic.ToneModulation" ToneModulation.tests : + prefix "Plain.ToneModulation" ToneModulationL.tests : + prefix "Generic.ToneModulation" ToneModulationG.tests : + prefix "Generic.Permutation" Permutation.tests : + prefix "Generic.Fourier" Fourier.tests : + prefix "Basic.NumberTheory" NumberTheory.tests : + prefix "Generic.FourierInteger" FourierInteger.tests : + prefix "Generic.Filter" FilterG.tests : + prefix "Causal.Analysis" AnalysisC.tests : + [] diff --git a/test/Test/Sound/Synthesizer/Basic/NumberTheory.hs b/test/Test/Sound/Synthesizer/Basic/NumberTheory.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Basic/NumberTheory.hs @@ -0,0 +1,119 @@ +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Sound.Synthesizer.Basic.NumberTheory (tests) where + +import Synthesizer.Basic.NumberTheory (Order(Order), ) +import qualified Synthesizer.Basic.NumberTheory as NT +import qualified Data.Set as Set + +import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) + +import qualified Algebra.Absolute as Absolute + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +newtype Cardinal a = Cardinal a + deriving (Show) + +instance (Absolute.C a, Arbitrary a) => Arbitrary (Cardinal a) where + arbitrary = fmap (Cardinal . abs) arbitrary + + +newtype Positive a = Positive a + deriving (Show) + +instance (Absolute.C a, Arbitrary a) => Arbitrary (Positive a) where + arbitrary = fmap (Positive . (1+) . abs) arbitrary + + +simple :: + (Testable t, Arbitrary (wrapper Integer), Show (wrapper Integer)) => + (wrapper Integer -> t) -> IO () +simple = quickCheck + +tests :: [(String, IO ())] +tests = + ("primitiveRootsOfUnity naive vs. power", + simple $ \(Cardinal m) order -> + NT.primitiveRootsOfUnityNaive m order + == + NT.primitiveRootsOfUnityPower m order) : + ("primitiveRootsOfUnity naive vs. fullorbit", + simple $ \(Cardinal m) order -> + NT.primitiveRootsOfUnityNaive m order + == + (Set.toAscList $ Set.fromList $ + NT.primitiveRootsOfUnityFullOrbit m order)) : + ("Carmichael theorem", + simple $ \(Positive a) (Positive b) -> + NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity (lcm a b)) + == + lcm + (NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity a)) + (NT.getOrder (NT.maximumOrderOfPrimitiveRootsOfUnity b))) : + ("maximumOrderOfPrimitiveRootsOfUnity naive vs. integer", + simple $ \(Positive m) -> + NT.maximumOrderOfPrimitiveRootsOfUnityNaive m + == + NT.maximumOrderOfPrimitiveRootsOfUnityInteger m) : + ("number of rootsOfUnityPower, lcm", + simple $ \(Positive m) ao@(Order a) bo@(Order b) -> + let g = length . NT.rootsOfUnityPower m + in g (Order $ lcm a b) == lcm (g ao) (g bo)) : + ("ringsWithPrimitiveRootsOfUnityAndUnits: minimal modulus", + quickCheck $ \order@(Order expo) -> + (head $ NT.ringsWithPrimitiveRootOfUnityAndUnit order) + == + (head $ NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive + [order] [expo])) : + ("combine two rings with primitive roots of certain orders", + quickCheck $ \m n -> + let r = lcm + (head (NT.ringsWithPrimitiveRootOfUnityAndUnit m)) + (head (NT.ringsWithPrimitiveRootOfUnityAndUnit n)) + in NT.hasPrimitiveRootOfUnityInteger r m + && + NT.hasPrimitiveRootOfUnityInteger r n) : + ("combine many rings with primitive roots of certain orders", + quickCheck $ \n0 ns0 -> + let ns = take 3 $ map (\n -> 1 + mod n 10) (n0:ns0) + order = NT.lcmMulti ns + in take 3 (NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive + (map Order ns) ns) + == + take 3 (NT.ringsWithPrimitiveRootsOfUnityAndUnitsNaive + [Order order] [order])) : +{- +Unfortunately rings with certain units cannot be combined +while maintaining these elements as units. + +Counterexample: + ringsWithPrimitiveRootOfUnityAndUnit 2 = 3:... + ringsWithPrimitiveRootOfUnityAndUnit 3 = 7:... + But in Z_{3·7} the number 3 is no unit. + + ("combine rings with certain units", + quickCheck $ \(Positive m) (Positive n) -> + let r = fromIntegral $ lcm + (head (NT.ringsWithPrimitiveRootOfUnityAndUnit m)) + (head (NT.ringsWithPrimitiveRootOfUnityAndUnit n)) + in PID.coprime r m && PID.coprime r n) : +-} + ("number of roots of unity lcm", + quickCheck $ \(Positive n) (Positive k) (Positive l) -> + let orders = NT.ordersOfRootsOfUnityInteger !! (n-1) + in lcm (orders!!(k-1)) (orders!!(l-1)) + == + orders !! (lcm k l - 1)) : + ("number of roots of unity vs. primitive roots", + quickCheck $ \(Positive n) (Positive k) -> + (sum $ map snd $ + filter (flip divides k . fst) $ + zip + [1..] + (NT.ordersOfPrimitiveRootsOfUnityInteger !! (n-1))) + == + NT.ordersOfRootsOfUnityInteger !! (n-1) !! (k-1)) : + [] diff --git a/test/Test/Sound/Synthesizer/Basic/ToneModulation.hs b/test/Test/Sound/Synthesizer/Basic/ToneModulation.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Basic/ToneModulation.hs @@ -0,0 +1,93 @@ +module Test.Sound.Synthesizer.Basic.ToneModulation where + +import qualified Synthesizer.Interpolation as Interpolation +import Synthesizer.Interpolation (margin, ) + +import qualified Synthesizer.Basic.Phase as Phase +import qualified Synthesizer.Basic.ToneModulation as ToneMod + +import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest + +import Test.QuickCheck (quickCheck, Property, (==>), Testable, ) +-- import Test.Utility + +import qualified Number.NonNegative as NonNeg + +import qualified Algebra.RealField as RealField +import qualified Algebra.Field as Field + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +untangleShapePhase :: (Field.C a, Eq a) => + Int -> a -> (a, a) -> Property +untangleShapePhase periodInt period c = + period /= zero ==> + ToneMod.untangleShapePhase periodInt period c == + ToneMod.untangleShapePhaseAnalytic periodInt period c + +flattenShapePhase :: (RealField.C a) => + Int -> a -> (a, Phase.T a) -> Property +flattenShapePhase periodInt period c = + period /= zero ==> + ToneMod.flattenShapePhase periodInt period c == + ToneMod.flattenShapePhaseAnalytic periodInt period c + + +-- * auxiliary quickCheck functions + +{- +Although that looks like a too small value, it is actually right, +because numberLeap counts intervals of size periodInt, not single elements. +So numberLeap=2 like in linear interpolation means 2*periodInt. +-} +minLength :: + Interpolation.T a v -> + Interpolation.T a v -> + Int -> NonNeg.Int -> Int +minLength ipLeap ipStep = + minLengthMargin (margin ipLeap) (margin ipStep) + +minLengthMargin :: + Interpolation.Margin -> + Interpolation.Margin -> + Int -> NonNeg.Int -> Int +minLengthMargin marginLeap marginStep periodInt ext = + ToneMod.interpolationNumber + marginLeap marginStep periodInt + + NonNeg.toNumber ext + + + +shapeLimits :: + Interpolation.T a v -> + Interpolation.T a v -> + Int -> Int -> (Int, Int) +shapeLimits ipLeap ipStep periodInt len = + ToneMod.shapeLimits + (margin ipLeap) (margin ipStep) + periodInt len + + + +testRationalLineIp :: Testable quickCheck => + (InterpolationTest.LinePreserving Rational Rational -> quickCheck) -> IO () +testRationalLineIp f = quickCheck f + +testRationalIp :: Testable quickCheck => + (InterpolationTest.T Rational Rational -> quickCheck) -> IO () +testRationalIp f = quickCheck f + + +tests :: [(String, IO ())] +tests = + ("untangleShapePhase", + quickCheck $ \periodInt period -> + untangleShapePhase periodInt (period :: Rational)) : + ("flattenShapePhase", + quickCheck $ \periodInt period -> + flattenShapePhase periodInt (period :: Rational)) : + [] diff --git a/test/Test/Sound/Synthesizer/Causal/Analysis.hs b/test/Test/Sound/Synthesizer/Causal/Analysis.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Causal/Analysis.hs @@ -0,0 +1,32 @@ +module Test.Sound.Synthesizer.Causal.Analysis (tests) where + +import qualified Synthesizer.Causal.Analysis as AnaC +import qualified Synthesizer.Causal.Process as Causal +import qualified Synthesizer.Plain.Analysis as Ana + +import Control.Arrow ((<<<), ) + +import qualified Data.List.Match as Match + +import Test.QuickCheck (quickCheck, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +tests :: [(String, IO ())] +tests = + ("deltaSigmaModulation", + quickCheck $ \xs -> + Match.take xs (Ana.deltaSigmaModulation xs) + == + Causal.apply AnaC.deltaSigmaModulation (xs::[Rational])) : + ("deltaSigmaModulationPositive", + quickCheck $ \threshold xs -> + Match.take xs (Ana.deltaSigmaModulationPositive threshold xs) + == + Causal.apply + (AnaC.deltaSigmaModulationPositive <<< + Causal.feedConstFst threshold) (xs::[Rational])) : + [] diff --git a/test/Test/Sound/Synthesizer/Generic/Cut.hs b/test/Test/Sound/Synthesizer/Generic/Cut.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/Cut.hs @@ -0,0 +1,104 @@ +module Test.Sound.Synthesizer.Generic.Cut (tests) where + +import qualified Synthesizer.Generic.Cut as CutG +import qualified Synthesizer.Generic.Signal as SigG + +import qualified Synthesizer.Storable.Signal as SigSt + +import qualified Synthesizer.ChunkySize.Signal as SigChunky +import qualified Synthesizer.ChunkySize as ChunkySize + +import qualified Data.StorableVector as SV +import qualified Data.StorableVector.Lazy.Pattern as SVP + +import qualified Synthesizer.State.Signal as SigS + +import qualified Data.EventList.Relative.BodyTime as EventList + +import qualified Number.NonNegative as NonNeg +import qualified Number.NonNegativeChunky as Chunky + +import qualified Numeric.NonNegative.Wrapper as NonNeg98 + +import Data.Tuple.HT (mapSnd, ) + +import Test.QuickCheck (quickCheck, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +dropMarginRemLength :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool +dropMarginRemLength nn nm xs = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + (k,ys) = CutG.dropMarginRem n m xs + in length xs - m == length ys - k + +dropMarginRemState :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool +dropMarginRemState nn nm xs = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + in CutG.dropMarginRem n m (SigS.fromList xs) + == + mapSnd SigS.fromList (CutG.dropMarginRem n m xs) + +dropMarginRemSV :: NonNeg.Int -> NonNeg.Int -> [Int] -> Bool +dropMarginRemSV nn nm xs = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + in CutG.dropMarginRem n m (SV.pack xs) + == + mapSnd SV.pack (CutG.dropMarginRem n m xs) + +dropMarginRemSVL :: NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> [Int] -> Bool +dropMarginRemSVL nn nm pat xs = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + in CutG.dropMarginRem n m + (CutG.take (CutG.length pat) xs) + == + mapSnd SigG.toList + (CutG.dropMarginRem n m + (SigChunky.fromState pat $ + SigG.toState xs :: SigSt.T Int)) + +dropMarginRemChunkySize :: + NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> Int -> Bool +dropMarginRemChunkySize nn nm pat x = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + in CutG.dropMarginRem n m pat + == + mapSnd + (ChunkySize.fromStorableVectorSize . SVP.length) + (CutG.dropMarginRem n m + (SVP.replicate (ChunkySize.toStorableVectorSize pat) x)) + +dropMarginRemPiecewise :: + NonNeg.Int -> NonNeg.Int -> ChunkySize.T -> Int -> Bool +dropMarginRemPiecewise nn nm pat x = + let n = NonNeg.toNumber nn + m = NonNeg.toNumber nm + in CutG.dropMarginRem n m pat + == + mapSnd + (Chunky.fromChunks . + map (\size -> SigG.LazySize $ NonNeg98.toNumber size) . + EventList.getTimes) + (CutG.dropMarginRem n m + (EventList.fromPairList $ map ((,) x) $ + map (\(SigG.LazySize size) -> NonNeg98.fromNumber size) $ + Chunky.toChunks pat)) + + +tests :: [(String, IO ())] +tests = + ("dropMarginRemLength", quickCheck dropMarginRemLength) : + ("dropMarginRemState", quickCheck dropMarginRemState) : + ("dropMarginRemSV", quickCheck dropMarginRemSV) : + ("dropMarginRemSVL", quickCheck dropMarginRemSVL) : + ("dropMarginRemChunkySize", quickCheck dropMarginRemChunkySize) : + ("dropMarginRemPiecewise", quickCheck dropMarginRemPiecewise) : + [] diff --git a/test/Test/Sound/Synthesizer/Generic/Filter.hs b/test/Test/Sound/Synthesizer/Generic/Filter.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/Filter.hs @@ -0,0 +1,64 @@ +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Sound.Synthesizer.Generic.Filter (tests) where + +import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG +import qualified Synthesizer.Generic.Cyclic as Cyclic +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.Generic.Cut as CutG +import qualified Synthesizer.Plain.Signal as Sig + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty + +import Test.QuickCheck (Testable, quickCheck, ) + +-- import qualified Algebra.Ring as Ring + +import qualified Algebra.Laws as Law + +import NumericPrelude.Numeric +import NumericPrelude.Base + + +simple :: + (Testable t) => + (Sig.T Int -> t) -> IO () +simple = quickCheck + +(=|=) :: + (Eq sig, CutG.Transform sig) => + sig -> sig -> Bool +x =|= y = + CutG.take 100 x == CutG.take 100 y + +tests :: [(String, IO ())] +tests = + ("identity", + simple $ Law.identity FiltNRG.generic $ SigG.singleton one) : + ("commutativity", + simple $ Law.commutative FiltNRG.generic) : + ("distributivity", + simple $ Law.leftDistributive FiltNRG.generic SigG.mix) : + ("karatsuba finite", + simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaFinite (*) x y) : + ("karatsuba finite-infinite", + simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaFiniteInfinite (*) x y) : + ("karatsuba infinite", + simple $ \x y -> FiltNRG.generic x y == FiltNRG.karatsubaInfinite (*) x y) : + ("karatsuba finite-infinite cycle", + simple $ \x yn -> + case NonEmpty.toInfiniteList yn of + y -> FiltNRG.generic x y =|= FiltNRG.karatsubaFiniteInfinite (*) x y) : + ("karatsuba infinite cycle", + simple $ \x yn -> + case NonEmpty.toInfiniteList yn of + y -> FiltNRG.generic x y =|= FiltNRG.karatsubaInfinite (*) x y) : + ("convolve triple", + quickCheck $ \x y -> + Cyclic.sumAndConvolveTriple x y == + Cyclic.sumAndConvolveTripleAlt x (y :: Cyclic.Triple Integer)) : + ("periodic summation", + simple $ \x y n -> + let periodic = Cyclic.fromSignal SigG.defaultLazySize (1 + abs n) + in Cyclic.convolve (periodic x) (periodic y) == + periodic (FiltNRG.generic x y)) : + [] diff --git a/test/Test/Sound/Synthesizer/Generic/Fourier.hs b/test/Test/Sound/Synthesizer/Generic/Fourier.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/Fourier.hs @@ -0,0 +1,151 @@ +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Sound.Synthesizer.Generic.Fourier (tests) where + +import qualified Synthesizer.Generic.Fourier as Fourier +import qualified Synthesizer.Generic.Cyclic as Cyclic +import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG +import qualified Synthesizer.Generic.Analysis as AnaG +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.Generic.Cut as CutG +import qualified Synthesizer.Storable.Signal as SigSt +import qualified Synthesizer.State.Signal as SigS + +import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) +import Test.Utility (approxEqualAbs, approxEqualComplexAbs, ) + +import qualified Number.Complex as Complex + +import qualified Algebra.Ring as Ring +import qualified Algebra.Additive as Additive + +import Control.Monad (liftM2, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +tolerance :: Double +tolerance = 1e-10 + +normalize :: + SigSt.T (Complex.T Double) -> SigSt.T (Complex.T Double) +normalize xs = + FiltNRG.amplifyVector + (recip $ max (0.1::Double) $ AnaG.volumeVectorMaximum xs) xs + +newtype Normed = Normed (SigSt.T (Complex.T Double)) + deriving (Show) + +instance Arbitrary Normed where + arbitrary = fmap (Normed . normalize) arbitrary + + +data Normed2 = + Normed2 + (SigSt.T (Complex.T Double)) + (SigSt.T (Complex.T Double)) + deriving (Show) + +instance Arbitrary Normed2 where + arbitrary = + liftM2 + (\x y -> + let len = min (CutG.length x) (CutG.length y) + in Normed2 + (normalize $ CutG.take len x) + (normalize $ CutG.take len y)) + arbitrary + arbitrary + + +-- could be moved to NumericPrelude +class Complex a where + conjugate :: a -> a + +instance (Additive.C a) => Complex (Complex.T a) where + conjugate = Complex.conjugate + +scalarProduct :: + (SigG.Read sig y, Ring.C y, Complex y) => + sig y -> sig y -> y +scalarProduct xs ys = + SigS.sum $ + SigS.zipWith (*) + (SigG.toState xs) + (SigS.map conjugate $ SigG.toState ys) + +(=~=) :: + SigSt.T (Complex.T Double) -> + SigSt.T (Complex.T Double) -> + Bool +(=~=) xs ys = + SigG.length xs == SigG.length ys && + (SigG.foldR (&&) True $ + SigG.zipWith (approxEqualComplexAbs tolerance) xs ys) + +simple :: + (Testable t) => + (SigSt.T (Complex.T Double) -> t) -> IO () +simple = quickCheck + +tests :: [(String, IO ())] +tests = + ("fourier inverse", + quickCheck $ \(Normed x) -> + x =~= + (FiltNRG.amplify (recip $ fromIntegral $ SigG.length x) $ + Fourier.transformBackward $ Fourier.transformForward x)) : + ("double fourier = reverse", + quickCheck $ \(Normed x) -> + x =~= + (Cyclic.reverse $ + FiltNRG.amplify (recip $ fromIntegral $ SigG.length x) $ + Fourier.transformForward $ + Fourier.transformForward x)) : + ("fourier of reverse", + quickCheck $ \(Normed x) -> + Cyclic.reverse (Fourier.transformForward x) =~= + Fourier.transformForward (Cyclic.reverse x)) : + ("fourier of conjugate", + quickCheck $ \(Normed x) -> + (SigG.map Complex.conjugate $ Fourier.transformForward x) + =~= + (Fourier.transformForward $ + SigG.map Complex.conjugate $ Cyclic.reverse x)) : + ("additivity", + quickCheck $ \(Normed2 x y) -> + SigG.mix (Fourier.transformForward x) (Fourier.transformForward y) + =~= + Fourier.transformForward (SigG.mix x y)) : + ("isometry", + simple $ \xs x0 -> + let x = normalize (SigG.cons x0 xs) + in approxEqualAbs tolerance + (AnaG.volumeVectorEuclideanSqr $ Fourier.transformForward x) + (fromIntegral (SigG.length x) * + AnaG.volumeVectorEuclideanSqr x)) : + ("unitarity", + quickCheck $ \(Normed2 x y) -> + approxEqualComplexAbs tolerance + (scalarProduct + (Fourier.transformForward x) (Fourier.transformForward y)) + (fromIntegral (SigG.length x) * scalarProduct x y)) : + ("convolution", + quickCheck $ \(Normed2 x y) -> + SigG.zipWith (*) + (Fourier.transformForward x) + (Fourier.transformForward y) + =~= + Fourier.transformForward (Cyclic.convolve x y)) : + ("convolution cyclic", + quickCheck $ \(Normed2 x y) -> + Fourier.convolveCyclic x y + =~= + Cyclic.convolve x y) : + ("convolution long", + quickCheck $ \(Normed x) (Normed y) -> + FiltNRG.karatsubaFinite (*) x y + =~= + Fourier.convolveWithWindow (Fourier.window x) y) : + [] diff --git a/test/Test/Sound/Synthesizer/Generic/FourierInteger.hs b/test/Test/Sound/Synthesizer/Generic/FourierInteger.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/FourierInteger.hs @@ -0,0 +1,178 @@ +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Sound.Synthesizer.Generic.FourierInteger (tests) where + +import qualified Synthesizer.Generic.Fourier as Fourier +import qualified Synthesizer.Generic.Cyclic as Cyclic +import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.Generic.Cut as CutG +import qualified Synthesizer.State.Signal as SigS +import qualified Synthesizer.Plain.Signal as Sig + +import Test.QuickCheck (Testable, Arbitrary, arbitrary, quickCheck, ) + +import qualified Synthesizer.Basic.NumberTheory as NT + +import qualified Number.ResidueClass.Check as RC +import Number.ResidueClass.Check ((/:), ) + +import qualified Algebra.ToInteger as ToInteger +import qualified Algebra.IntegralDomain as Integral +import qualified Algebra.Ring as Ring + +import Control.Monad (liftM2, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +newtype Modulus a = Modulus a + deriving (Show) + +instance Ring.C a => Arbitrary (Modulus a) where + arbitrary = fmap (Modulus . (2+) . fromInteger) arbitrary + + +data ModularSignal = + ModularSignal (Modulus Integer) (Sig.T (RC.T Integer)) + deriving (Show) + +instance Arbitrary ModularSignal where + arbitrary = + fmap (uncurry ModularSignal . signal) arbitrary + + +data ModularSignal2 = + ModularSignal2 + (Modulus Integer) (Sig.T (RC.T Integer)) (Sig.T (RC.T Integer)) + deriving (Show) + +instance Arbitrary ModularSignal2 where + arbitrary = + liftM2 + (\x y -> + let len = min (CutG.length x) (CutG.length y) + m = NT.fastFourierRing len + in ModularSignal2 + (Modulus m) + (fmap (/: m) $ CutG.take len x) + (fmap (/: m) $ CutG.take len y)) + arbitrary + arbitrary + +scalarProduct :: + Modulus Integer -> + Sig.T (RC.T Integer) -> Sig.T (RC.T Integer) -> + RC.T Integer +scalarProduct (Modulus m) xs ys = + SigS.foldL (+) (RC.zero m) $ + SigS.zipWith (*) + (SigG.toState xs) + (SigG.toState ys) + +{- +signal :: + Integral.C a => + Modulus a -> Sig.T a -> Sig.T (RC.T a) +signal (Modulus a) = fmap (/: a) +-} + +signal :: + Sig.T Integer -> (Modulus Integer, Sig.T (RC.T Integer)) +signal xs = + let m = NT.fastFourierRing $ length xs + in (Modulus m, fmap (/: m) xs) + +modular :: + (Integral.C a, ToInteger.C b) => + Modulus a -> b -> RC.T a +modular (Modulus m) = + RC.fromRepresentative m . fromIntegral + + +simple :: + (Testable t) => + (Sig.T Integer -> t) -> IO () +simple = quickCheck + +tests :: [(String, IO ())] +tests = + ("fourier inverse", + quickCheck $ \(ModularSignal m x) -> + (Fourier.transformBackward $ Fourier.transformForward x) + == + FiltNRG.amplify (modular m $ length x) x) : + ("double fourier = reverse", + quickCheck $ \(ModularSignal m x) -> + (Cyclic.reverse $ + Fourier.transformForward $ + Fourier.transformForward x) + == + FiltNRG.amplify (modular m $ length x) x) : + ("fourier of reverse", + quickCheck $ \(ModularSignal _m x) -> + Cyclic.reverse (Fourier.transformForward x) == + Fourier.transformForward (Cyclic.reverse x)) : + ("homogenity", + quickCheck $ \(ModularSignal m x) y -> + (FiltNRG.amplify (modular m (y::Integer)) $ + Fourier.transformForward x) + == + (Fourier.transformForward $ + FiltNRG.amplify (modular m y) x)) : + ("additivity", + quickCheck $ \(ModularSignal2 _m x y) -> + SigG.mix (Fourier.transformForward x) (Fourier.transformForward y) + == + Fourier.transformForward (SigG.mix x y)) : +{- + ("isometry", + simple $ \xs x0 -> + let (m,x) = signal (SigG.cons x0 xs) + in (AnaG.volumeVectorEuclideanSqr $ Fourier.transformForward x) + == + (modular m (SigG.length x) * + AnaG.volumeVectorEuclideanSqr x)) : +-} + ("unitarity", + quickCheck $ \(ModularSignal2 m x y) -> + {- + since there is no equivalent of a complex conjugate + we have to take the scalar product with the backwards transform. + -} + scalarProduct m + (Fourier.transformForward x) (Fourier.transformBackward y) + == + modular m (length x) * scalarProduct m x y) : + ("convolution", + quickCheck $ \(ModularSignal2 _m x y) -> + SigG.zipWith (*) + (Fourier.transformForward x) + (Fourier.transformForward y) + == + Fourier.transformForward (Cyclic.convolve x y)) : + ("convolution cyclic", + quickCheck $ \(ModularSignal2 _m x y) -> + Fourier.convolveCyclic x y + == + Cyclic.convolve x y) : + ("convolution long", + simple $ \x0 y0 -> + let m = Modulus $ NT.fastFourierRing $ + 2 * (NT.ceilingPowerOfTwo $ length x0) + x = fmap (modular m) x0 + y = fmap (modular m) y0 + in fmap (modular m) (FiltNRG.karatsubaFinite (*) x0 y0) + == + Fourier.convolveWithWindow (Fourier.window x) y) : + ("convolution long modular", + simple $ \x0 y0 -> + let m = Modulus $ NT.fastFourierRing $ + 2 * (NT.ceilingPowerOfTwo $ length x0) + x = fmap (modular m) x0 + y = fmap (modular m) (y0 :: Sig.T Integer) + in FiltNRG.karatsubaFinite (*) x y + == + Fourier.convolveWithWindow (Fourier.window x) y) : + [] diff --git a/test/Test/Sound/Synthesizer/Generic/Permutation.hs b/test/Test/Sound/Synthesizer/Generic/Permutation.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/Permutation.hs @@ -0,0 +1,45 @@ +{- +wish list: + - custom Permutation type with Arbitrary instance +-} +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Sound.Synthesizer.Generic.Permutation (tests) where + +import qualified Synthesizer.Generic.Permutation as Permutation + +import Test.QuickCheck (quickCheck, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +tests :: [(String, IO ())] +tests = + ("inverse transposition", + quickCheck $ \n0 m0 -> + let n = mod n0 100 + m = mod m0 100 + in Permutation.inverse (Permutation.transposition n m) + == + Permutation.transposition m n) : + ("inverse skewGrid", + quickCheck $ \n0 m0 -> + let g = gcd n0 m0 + (n,m) = if g==0 then (0,0) else (abs (div n0 g), abs (div m0 g)) + in Permutation.inverse (Permutation.skewGrid n m) + == + Permutation.skewGridInv n m) : + ("inverse skewGridCRT", + quickCheck $ \n0 m0 -> + let g = gcd n0 m0 + (n,m) = if g==0 then (0,0) else (abs (div n0 g), abs (div m0 g)) + in Permutation.inverse (Permutation.skewGridCRT n m) + == + Permutation.skewGridCRTInv n m) : + {- + reverse (multiplicative (generator n) n) + == + multiplicative (recip $ generator n) n + -} + [] diff --git a/test/Test/Sound/Synthesizer/Generic/ToneModulation.hs b/test/Test/Sound/Synthesizer/Generic/ToneModulation.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Generic/ToneModulation.hs @@ -0,0 +1,304 @@ +module Test.Sound.Synthesizer.Generic.ToneModulation (tests) where + +import Test.Sound.Synthesizer.Basic.ToneModulation ( + minLength, + minLengthMargin, +-- shapeLimits, +-- testRationalLineIp, + testRationalIp, + ) + +import qualified Synthesizer.Causal.ToneModulation as ToneModC +import qualified Synthesizer.Generic.Wave as WaveG + +import qualified Synthesizer.Plain.Signal as Sig +import qualified Synthesizer.Plain.Oscillator as Osci +import qualified Synthesizer.Plain.Interpolation as Interpolation +import qualified Synthesizer.Plain.ToneModulation as ToneModL +import qualified Synthesizer.Plain.Wave as WaveL +import Synthesizer.Interpolation (marginNumber, ) + +import qualified Synthesizer.Causal.Oscillator as OsciC +import qualified Synthesizer.Causal.Process as Causal + +import qualified Synthesizer.State.Signal as SigS + +import qualified Synthesizer.Basic.Wave as Wave +import qualified Synthesizer.Basic.Phase as Phase + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty +import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest + +import Test.QuickCheck (quickCheck, Property, (==>), ) +import Test.Utility (ArbChar, ) +-- import Debug.Trace (trace, ) + +import qualified Number.NonNegative as NonNeg + +import qualified Algebra.RealField as RealField + + +import Data.List.HT (viewL, takeWhileJust, ) +import Data.Tuple.HT (mapSnd, ) +import qualified Data.List as List + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +limitMinRelativeValues :: + Int -> Int -> [NonNeg.Int] -> Bool +limitMinRelativeValues xMin x0 xsnn = + let xs = map NonNeg.toNumber xsnn + (y0,limiter) = ToneModC.limitMinRelativeValues xMin x0 + in (y0, Causal.apply limiter xs) == + ToneModL.limitMinRelativeValues xMin x0 xs + +integrateFractional :: (RealField.C t) => + NonNeg.T t -> t -> Phase.T t -> [NonNeg.T t] -> [t] -> Property +integrateFractional + periodNN shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + (c0, coordinator) = + ToneModC.integrateFractional + period (shape0, phase) + coords = + ToneModL.integrateFractional + period (shape0, shapes) (phase, freqs) + in period /= zero ==> + c0 : Causal.apply coordinator (zip shapes freqs) == + coords + +-- oscillatorCellSize :: (Show t, Show v, RealField.C t, Eq v) => +oscillatorCellSize :: (RealField.C t, Eq v) => + Interpolation.Margin -> + Interpolation.Margin -> + NonNeg.Int -> NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> + Property +oscillatorCellSize + marginLeap marginStep periodIntNN periodNN ext + ixs shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + len = minLengthMargin marginLeap marginStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledTone = + ToneModC.oscillatorCells + marginLeap marginStep periodInt period tone + (shape0, Phase.fromRepresentative phase) + `Causal.apply` + zip shapes freqs + in period /= zero && + marginNumber marginLeap > zero && + marginNumber marginStep > zero ==> + all + ((\cell -> + Sig.lengthAtLeast (marginNumber marginLeap) cell && + all (Sig.lengthAtLeast (marginNumber marginStep)) + (take (marginNumber marginLeap) cell)) + . SigS.toList . snd) + resampledTone + +oscillatorSuffixes :: (RealField.C t, Eq v) => + Interpolation.Margin -> + Interpolation.Margin -> + NonNeg.Int -> NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> + Property +oscillatorSuffixes + marginLeap marginStep periodIntNN periodNN ext + ixs shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + len = minLengthMargin marginLeap marginStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledToneA = + init $ + map (\(sp,cell) -> + (sp, takeWhileJust . map (fmap fst . viewL) $ cell)) $ + ToneModL.oscillatorSuffixes + marginLeap marginStep periodInt period tone + (shape0, shapes) (Phase.fromRepresentative phase, freqs) + resampledToneB = + ToneModC.oscillatorSuffixes + marginLeap marginStep periodInt period tone + (shape0, Phase.fromRepresentative phase) + `Causal.apply` + zip shapes freqs + in period /= zero && + periodInt /= zero && + marginNumber marginLeap > zero && + marginNumber marginStep > zero ==> + resampledToneA == resampledToneB + +oscillatorCells :: (RealField.C t, Eq v) => + Interpolation.Margin -> + Interpolation.Margin -> + NonNeg.Int -> NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> + Property +oscillatorCells + marginLeap marginStep periodIntNN periodNN ext + ixs shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + len = minLengthMargin marginLeap marginStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledToneA = + init $ map (mapSnd List.transpose) $ + ToneModL.oscillatorCells + marginLeap marginStep periodInt period tone + (shape0, shapes) (Phase.fromRepresentative phase, freqs) + resampledToneB = + map (mapSnd SigS.toList) $ + ToneModC.oscillatorCells + marginLeap marginStep periodInt period tone + (shape0, Phase.fromRepresentative phase) + `Causal.apply` + zip shapes freqs + in period /= zero && + periodInt /= zero && + marginNumber marginLeap > zero && + marginNumber marginStep > zero ==> + resampledToneA == resampledToneB +{- +Margin {marginNumber = 1, marginOffset = 2} +Margin {marginNumber = 5, marginOffset = 5} +3 % 4 +0 +('\DEL',['~','~','"']) +-2 % 5 +2 % 5 +[2 % 1,3 % 4] +[-5 % 2,-1 % 2] +-} + +{- | +'WaveL.sampledTone' and 'WaveG.sampledTone' +do not only differ in the signal types they process, +but also in the way they order the signal values. +The cells for 'WaveL.sampledTone' are transposed +with respect to 'WaveG.sampledTone'. +-} +sampledTone :: (RealField.C a, Eq v) => + InterpolationTest.T a v -> + InterpolationTest.T a v -> + NonNeg.T a -> NonNeg.Int -> NonEmpty.T v -> + a -> Phase.T a -> Property +sampledTone = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs shape phase -> + let period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + in period /= zero ==> + WaveG.sampledTone ipLeap ipStep period tone shape `Wave.apply` phase == + WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` phase + + + +shapeFreqModFromSampledTone :: (RealField.C t, Eq v) => + InterpolationTest.T t v -> + InterpolationTest.T t v -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> Phase.T t -> [NonNeg.T t] -> [t] -> + Property +shapeFreqModFromSampledTone = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs shape0 phase shapesNN freqs -> + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledToneA = + init $ + Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone + shape0 (Phase.toRepresentative phase) shapes freqs + resampledToneB = + OsciC.shapeFreqModFromSampledTone + ipLeap ipStep period tone shape0 phase + `Causal.apply` + zip shapes freqs + in period /= zero ==> + resampledToneA == resampledToneB + + +{- +We have a problem here with the phase distortion signal, +because frequency and shape modulation control signals +are delayed by one element with respect to the phase distortion. +How can we cope with different lengths of the control signals, +without padding the phase control with zeros? +This one did not work: + phaseDistorts = pd:pds + resampledToneA = + Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone + shape0 (Phase.toRepresentative phase) shapes (init phaseDistorts) freqs +-} +shapePhaseFreqModFromSampledTone :: (RealField.C t, Eq v) => + InterpolationTest.T t v -> + InterpolationTest.T t v -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> Phase.T t -> [NonNeg.T t] -> (t,[t]) -> [t] -> + Property +shapePhaseFreqModFromSampledTone = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs shape0 phase shapesNN (pd,pds) freqs -> + let period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + shapes = map NonNeg.toNumber shapesNN + phaseDistorts = pd:pds ++ repeat zero + resampledToneA = + init $ + Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone + shape0 (Phase.toRepresentative phase) shapes phaseDistorts freqs + resampledToneB = + OsciC.shapePhaseFreqModFromSampledTone + ipLeap ipStep period tone shape0 phase + `Causal.apply` + zip3 shapes phaseDistorts freqs + in period /= zero ==> + resampledToneA == resampledToneB + + + +tests :: [(String, IO ())] +tests = + ("limitMinRelativeValues", quickCheck limitMinRelativeValues) : + ("integrateFractional", + quickCheck (\period -> integrateFractional (period :: NonNeg.Rational))) : + ("oscillatorCellSize", + quickCheck (\ml ms periodInt period ext ixs -> + oscillatorCellSize ml ms periodInt (period :: NonNeg.Rational) + ext (ixs :: NonEmpty.T ArbChar))) : + ("oscillatorSuffixes", + quickCheck (\ml ms periodInt period ext ixs -> + oscillatorSuffixes ml ms periodInt (period :: NonNeg.Rational) + ext (ixs :: NonEmpty.T ArbChar))) : + ("oscillatorCells", + quickCheck (\ml ms periodInt period ext ixs -> + oscillatorCells ml ms periodInt (period :: NonNeg.Rational) + ext (ixs :: NonEmpty.T ArbChar))) : + ("sampledTone", + testRationalIp sampledTone) : + ("shapeFreqModFromSampledTone", + testRationalIp shapeFreqModFromSampledTone) : + ("shapePhaseFreqModFromSampledTone", + testRationalIp shapePhaseFreqModFromSampledTone) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Analysis.hs b/test/Test/Sound/Synthesizer/Plain/Analysis.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Analysis.hs @@ -0,0 +1,160 @@ +module Test.Sound.Synthesizer.Plain.Analysis (tests) where + +import qualified Synthesizer.Plain.Analysis as Analysis + +import qualified Algebra.Algebraic as Algebraic +import qualified Algebra.RealField as RealField +import qualified Algebra.Field as Field +import qualified Algebra.RealRing as RealRing + +import qualified Algebra.NormedSpace.Maximum as NormedMax +import qualified Algebra.NormedSpace.Euclidean as NormedEuc +import qualified Algebra.NormedSpace.Sum as NormedSum + +import qualified MathObj.LaurentPolynomial as LPoly + +import qualified Data.NonEmpty as NonEmpty +import Data.List (genericLength) + +import Test.QuickCheck (quickCheck, Property, (==>)) +import Test.Utility (approxEqual) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +volumeVectorMaximum :: (NormedMax.C y y, RealRing.C y) => [y] -> Bool +volumeVectorMaximum xs = + Analysis.volumeVectorMaximum xs == Analysis.volumeMaximum xs + +volumeVectorEuclidean :: + (NormedEuc.C y y, Algebraic.C y, Eq y) => + NonEmpty.T [] y -> Bool +volumeVectorEuclidean xs = + let ys = NonEmpty.flatten xs + in Analysis.volumeVectorEuclidean ys == Analysis.volumeEuclidean ys + +volumeVectorEuclideanSqr :: + (NormedEuc.Sqr y y, Field.C y, Eq y) => + NonEmpty.T [] y -> Bool +volumeVectorEuclideanSqr xs = + let ys = NonEmpty.flatten xs + in Analysis.volumeVectorEuclideanSqr ys == Analysis.volumeEuclideanSqr ys + +volumeVectorSum :: + (NormedSum.C y y, RealField.C y) => + NonEmpty.T [] y -> Bool +volumeVectorSum xs = + let ys = NonEmpty.flatten xs + in Analysis.volumeVectorSum ys == Analysis.volumeSum ys + + + +bounds :: Ord a => NonEmpty.T [] a -> Bool +bounds xs = + Analysis.bounds xs == (NonEmpty.minimum xs, NonEmpty.maximum xs) + + +spread :: RealField.C a => (a,a) -> Bool +spread b = + sum (map snd (Analysis.spread b)) == one + + +histogramDiscrete :: NonEmpty.T [] Int -> Bool +histogramDiscrete xs = + Analysis.histogramDiscreteArray xs == + Analysis.histogramDiscreteIntMap xs + +withEmptyHistogram :: + (NonEmpty.T [] y -> (Int, [y])) -> + [y] -> (Int, [y]) +withEmptyHistogram f = + maybe (error "no bounds", []) f . NonEmpty.fetch + +histogramDiscreteLength :: [Int] -> Bool +histogramDiscreteLength xs = + sum (snd (withEmptyHistogram Analysis.histogramDiscreteIntMap xs)) + == + length xs + +histogramDiscreteConcat :: [Int] -> [Int] -> Bool +histogramDiscreteConcat xs ys = + let xHist = withEmptyHistogram Analysis.histogramDiscreteIntMap xs + yHist = withEmptyHistogram Analysis.histogramDiscreteIntMap ys + xyHist0 = + LPoly.add + (uncurry LPoly.Cons xHist) + (uncurry LPoly.Cons yHist) + xyHist1 = + uncurry LPoly.Cons + (withEmptyHistogram Analysis.histogramDiscreteIntMap (xs++ys)) + in if null (LPoly.coeffs xyHist0) + then LPoly.coeffs xyHist0 == LPoly.coeffs xyHist1 + else xyHist0 == xyHist1 + + +histogramLinear :: NonEmpty.T [] Int -> Bool +histogramLinear xs = + let ys = fmap fromIntegral xs :: NonEmpty.T [] Double + in Analysis.histogramLinearArray ys == + Analysis.histogramLinearIntMap ys + + +histogramLinearLength :: NonEmpty.T [] Int -> Bool +histogramLinearLength xs = + let ys = fmap fromIntegral xs :: NonEmpty.T [] Double + in approxEqual 1e-10 + (genericLength $ NonEmpty.tail ys) + (sum (snd (Analysis.histogramLinearIntMap ys))) +{- +With eps = 1e-15 + +Falsifiable, after 83 tests: +-20 +[32,-41,11,-25,-17,-27,32,-36,7,-36,38] + +Falsifiable, after 78 tests: +10 +[-35,-28,-28,-24,-4,-29,-14,-29,-20,7,33,-2,-14,-4,7,-40,-5,-12] +-} + + + +centroid :: (Field.C a, Eq a) => [a] -> Property +centroid xs = + sum xs /= zero ==> + Analysis.centroid xs == Analysis.centroidAlt xs +-- Test.QuickCheck.quickCheck (\xs -> sum xs /= 0 Test.QuickCheck.==> propCentroid (xs::[Rational])) + +histogramDCOffset :: NonEmpty.T (NonEmpty.T []) Int -> Property +histogramDCOffset xs = + let x1 = NonEmpty.flatten xs + x = NonEmpty.flatten x1 + (offset, hist) = Analysis.histogramDiscreteArray x1 + in sum x /= 0 ==> + fromIntegral offset + Analysis.centroid (map fromIntegral hist) == + (Analysis.directCurrentOffset (map fromIntegral x) :: Rational) + + +small :: (Functor f) => f Int -> f Int +small = fmap (flip mod 1000) + + +tests :: [(String, IO ())] +tests = + ("volumeVectorMaximum", quickCheck (volumeVectorMaximum :: [Rational] -> Bool)) : + -- quickCheck may fail due to rounding errors, but so far the computation is exactly the same + ("volumeVectorEuclidean", quickCheck (volumeVectorEuclidean :: NonEmpty.T [] Double -> Bool)) : + ("volumeVectorEuclideanSqr", quickCheck (volumeVectorEuclideanSqr :: NonEmpty.T [] Rational -> Bool)) : + ("volumeVectorSum", quickCheck (volumeVectorSum :: NonEmpty.T [] Rational -> Bool)) : + ("bounds", quickCheck (bounds :: NonEmpty.T [] Rational -> Bool)) : + ("spread", quickCheck (spread :: (Rational,Rational) -> Bool)) : + ("histogramDiscrete", quickCheck (histogramDiscrete . small)) : + ("histogramDiscreteLength", quickCheck (histogramDiscreteLength . small)) : + ("histogramDiscreteConcat", quickCheck (\x y -> histogramDiscreteConcat (small x) (small y))) : + ("histogramLinear", quickCheck (histogramLinear . small)) : + ("histogramLinearLength", quickCheck (histogramLinearLength . small)) : + ("centroid", quickCheck (centroid :: [Rational] -> Property)) : + ("histogramDCOffset", quickCheck (histogramDCOffset . small)) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Control.hs b/test/Test/Sound/Synthesizer/Plain/Control.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Control.hs @@ -0,0 +1,112 @@ +module Test.Sound.Synthesizer.Plain.Control (tests) where + +import qualified Synthesizer.Plain.Control as Control + +import Test.QuickCheck (quickCheck, Property, (==>)) +import Test.Utility (equalList, approxEqualListAbs, approxEqualListRel, ) + +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive + +import Data.List (transpose) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +linearRing :: Int -> Int -> Bool +linearRing d y0 = +-- Control.linear d y0 == Control.linearMultiscale d y0 + all equalList $ take 100 $ transpose $ + Control.linear d y0 : + Control.linearMultiscale d y0 : + Control.linearStable d y0 : + [] + +{- +*Synthesizer.Plain.Control> propLinearApprox (-2/3) 2 +False + +Need a different definition of approximate equality. +-} +linearApprox :: Double -> Double -> Bool +linearApprox d y0 = + all (approxEqualListAbs (1e-10 * max (abs d) (abs y0))) $ + take 100 $ transpose $ + Control.linear d y0 : + Control.linearMean d y0 : + Control.linearMultiscale d y0 : + Control.linearStable d y0 : + [] + +linearExact :: Rational -> Rational -> Bool +linearExact d y0 = + all equalList $ take 100 $ transpose $ + Control.linear d y0 : + Control.linearMean d y0 : + Control.linearMultiscale d y0 : + Control.linearStable d y0 : + [] + +{- +Plain.Control.exponential: Falsifiable, after 88 tests: +-8.333333333333326e-2 +3.375 + +Plain.Control.exponential: Falsifiable, after 69 tests: +9.090909090909083e-2 +-10.0 + +Plain.Control.exponential: Falsifiable, after 73 tests: +-0.125 +-1.1428571428571428 + +Plain.Control.exponential2: Falsifiable, after 33 tests: +-7.692307692307687e-2 +9.5 +-} +exponential :: Double -> Double -> Bool +exponential time y0 = + all (approxEqualListRel (1e-10)) $ take 100 $ transpose $ + Control.exponential time y0 : + Control.exponentialMultiscale time y0 : + Control.exponentialStable time y0 : + [] + +exponential2 :: Double -> Double -> Bool +exponential2 time y0 = + all (approxEqualListRel (1e-10)) $ take 100 $ transpose $ + Control.exponential2 time y0 : + Control.exponential2Multiscale time y0 : + Control.exponential2Stable time y0 : + [] + +cosine :: Double -> Double -> Property +cosine t0 t1 = t0/=t1 ==> + all (approxEqualListAbs (1e-10)) $ + take 100 $ transpose $ + Control.cosine t0 t1 : + Control.cosineMultiscale t0 t1 : + Control.cosineStable t0 t1 : + [] + + +cubic :: (Rational, (Rational, Rational)) -> + (Rational, (Rational, Rational)) -> Property +cubic node0 node1 = fst node0 /= fst node1 ==> + take 100 (Control.cubicHermite node0 node1) == + take 100 (Control.cubicHermiteStable node0 node1) + + + +tests :: [(String, IO ())] +tests = + ("linearRing", quickCheck linearRing) : + ("linearApprox", quickCheck linearApprox) : + ("linearExact", quickCheck linearExact) : + ("exponential", quickCheck exponential) : + ("exponential2", quickCheck exponential2) : + ("cosine", quickCheck cosine) : + ("cubic", quickCheck cubic) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Filter.hs b/test/Test/Sound/Synthesizer/Plain/Filter.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Filter.hs @@ -0,0 +1,199 @@ +module Test.Sound.Synthesizer.Plain.Filter (tests) where + +import qualified Synthesizer.Plain.Filter.Recursive.MovingAverage as MA +import qualified Synthesizer.Plain.Filter.NonRecursive as FiltNR +import qualified Synthesizer.Plain.Signal as Sig +import qualified Synthesizer.Generic.Filter.NonRecursive as FiltNRG +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.Storable.Filter.NonRecursive as FiltNRSt +import qualified Synthesizer.Storable.Signal as SigSt +import qualified Synthesizer.Causal.Filter.NonRecursive as FiltNRC +import qualified Synthesizer.Causal.Process as Causal +import qualified Synthesizer.Frame.Stereo as Stereo + +import qualified Data.StorableVector.Lazy.Pattern as VP + +import Foreign.Storable.Tuple () + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty + +import Test.QuickCheck (quickCheck, {- Property, (==>) -}) +import Test.Utility (equalList, ArbChar, ) + +-- import qualified Algebra.Module as Module +-- import qualified Algebra.RealField as RealField +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive + +import qualified Number.GaloisField2p32m5 as GF +import qualified Number.NonNegative as NonNeg + +import qualified Numeric.NonNegative.Chunky as Chunky + +import qualified Data.List as List +import Data.Tuple.HT (mapPair, ) + +-- import Debug.Trace (trace, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +sums :: NonNeg.Int -> Rational -> Sig.T Rational -> Bool +sums nn x0 xs0 = + let n = min (length xs) (1 + NonNeg.toNumber nn) + xs = x0:xs0 + naive = FiltNR.sums n xs + pyramid = FiltNR.sumsPyramid n xs + rec = drop (n-1) $ MA.sumsStaticInt n xs + in -- this checks only for equal prefixes and can easily go wrong, + -- if one list is empty + and $ zipWith3 (\x y z -> x==y && y==z) naive rec pyramid + -- equalList $ naive : pyramid : rec : [] + +sumRange :: NonNeg.Int -> (NonNeg.Int, NonNeg.Int) -> Sig.T Int -> Bool +sumRange nheight (nl,nr) xs = + let wrap n = mod (NonNeg.toNumber n) (length xs + 1) + height = 1 + NonNeg.toNumber nheight + rng = (wrap nl, wrap nr) + pyr = take height (FiltNR.pyramid xs) + pyrSt = + FiltNRSt.pyramid (+) height + (SigSt.fromList SigSt.defaultChunkSize xs) + in equalList $ + FiltNR.sumRange xs rng : + FiltNR.sumRangeFromPyramid pyr rng : + FiltNR.sumRangeFromPyramidRec pyr rng : + FiltNR.sumRangeFromPyramidFoldr pyr rng : + FiltNRG.sumRangeFromPyramid pyrSt rng : + FiltNRG.sumRangeFromPyramidFoldr pyrSt rng : + FiltNRG.sumRangeFromPyramidReverse pyrSt rng : + [] + +getRange :: (NonNeg.Int, NonNeg.Int) -> NonEmpty.T (NonEmpty.T ArbChar) -> Bool +getRange (nl,nr) pyr0 = + let l = NonNeg.toNumber nl + r = NonNeg.toNumber nr + rng = if l<=r then (l,r) else (r,l) + pyr = map NonEmpty.toInfiniteList $ NonEmpty.toList pyr0 + in equalList $ + FiltNR.getRangeFromPyramid pyr rng : + FiltNRG.consumeRangeFromPyramid (:) [] pyr rng : + [] + +sumsPosModulated :: + NonNeg.Int -> Sig.T (NonNeg.Int,NonNeg.Int) -> NonEmpty.T Int -> Bool +sumsPosModulated nheight nctrl xsc = + let ctrl = map (mapPair (NonNeg.toNumber, NonNeg.toNumber)) nctrl + xs = NonEmpty.toInfiniteList xsc + height = min 10 $ NonNeg.toNumber nheight + in -- trace (show (height, ctrl, xsc)) $ + equalList $ + FiltNR.sumsPosModulated ctrl xs : + FiltNR.sumsPosModulatedPyramid height ctrl xs : + FiltNRG.sumsPosModulatedPyramid height ctrl xs : + SigSt.toList + (FiltNRG.sumsPosModulatedPyramid + height + (SigSt.fromList SigSt.defaultChunkSize ctrl) + (SigSt.fromList SigSt.defaultChunkSize xs)) : + SigSt.toList + (FiltNRSt.sumsPosModulatedPyramid + height + (SigSt.fromList SigSt.defaultChunkSize ctrl) + (SigSt.fromList SigSt.defaultChunkSize xs)) : + Causal.apply + (FiltNRC.sumsPosModulatedFromPyramid $ + FiltNRSt.pyramid (+) height $ + SigSt.fromList SigSt.defaultChunkSize xs) + ctrl : + [] + +minPosModulated :: + NonNeg.Int -> Sig.T (NonNeg.Int,NonNeg.Int) -> NonEmpty.T Int -> Bool +minPosModulated nheight nctrl xsc = + let ctrl = + map (\(nl,nr) -> + if nl==nr + then (NonNeg.toNumber nl, NonNeg.toNumber nr+1) + else (NonNeg.toNumber nl, NonNeg.toNumber nr)) + nctrl + xs = NonEmpty.toInfiniteList xsc + height = min 10 $ NonNeg.toNumber nheight + in -- trace (show (height, ctrl, xsc)) $ + equalList $ + zipWith FiltNR.minRange (List.tails xs) ctrl : + SigSt.toList + (FiltNRSt.accumulateBinPosModulatedPyramid min height + (SigSt.fromList SigSt.defaultChunkSize ctrl) + (SigSt.fromList SigSt.defaultChunkSize xs)) : + [] + +downSample2 :: + [Int] -> (Int, Sig.T Int) -> Bool +downSample2 lazySize xsc = + let len = Chunky.fromChunks $ map (VP.chunkSize . succ . abs) lazySize + xs = VP.pack len $ cycle $ uncurry (:) xsc + in equalList $ + FiltNRG.downsample2 SigG.defaultLazySize xs : + FiltNRSt.downsample2 xs : + [] + +sumsDownSample2 :: + [Int] -> (Int, Sig.T Int) -> Bool +sumsDownSample2 lazySize xsc = + let len = Chunky.fromChunks $ map (VP.chunkSize . succ . abs) lazySize + xs = VP.pack len $ cycle $ uncurry (:) xsc + in equalList $ + FiltNRG.sumsDownsample2 SigG.defaultLazySize xs : + FiltNRSt.sumsDownsample2 xs : + FiltNRSt.sumsDownsample2Alt xs : + [] + +{- +sumsDownSample2 :: + [VP.ChunkSize] -> (Int, Sig.T Int) -> Bool +sumsDownSample2 lazySize xsc = + let len = Chunky.fromChunks $ filter (0/=) lazySize + xs = VP.pack len $ cycle $ uncurry (:) xsc + in equalList $ + FiltNRG.sumsDownsample2 SigG.defaultLazySize xs : + FiltNRSt.sumsDownsample2 xs : + FiltNRSt.sumsDownsample2Alt xs : + [] +-} + +movingAverageModulatedPyramid :: + NonNeg.Int -> Sig.T NonNeg.Int -> + (Stereo.T GF.T, Sig.T (Stereo.T GF.T)) -> Bool +movingAverageModulatedPyramid nheight nctrl xsc = + let ctrl = map NonNeg.toNumber nctrl + xs = uncurry (:) xsc + pack ys = SigSt.fromList SigSt.defaultChunkSize ys + maxC = maximum ctrl + height = min 10 $ NonNeg.toNumber nheight + onegf :: GF.T + onegf = one + in -- trace (show (height, ctrl, xsc)) $ + equalList $ + pack (FiltNR.movingAverageModulatedPyramid onegf + height maxC ctrl (cycle xs)) : + FiltNRG.movingAverageModulatedPyramid onegf + height maxC (pack ctrl) (SigG.cycle $ pack xs) : + FiltNRSt.movingAverageModulatedPyramid onegf + height maxC (pack ctrl) (SigG.cycle $ pack xs) : + [] + + +tests :: [(String, IO ())] +tests = + ("sums", quickCheck sums) : + ("sumRange", quickCheck sumRange) : + ("getRange", quickCheck getRange) : + ("sumsPosModulated", quickCheck sumsPosModulated) : + ("minPosModulated", quickCheck minPosModulated): + ("downSample2", quickCheck downSample2) : + ("sumsDownSample2", quickCheck sumsDownSample2) : + ("movingAverageModulatedPyramid", quickCheck movingAverageModulatedPyramid) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs b/test/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Filter/Allpass.hs @@ -0,0 +1,56 @@ +module Test.Sound.Synthesizer.Plain.Filter.Allpass (tests) where + +import qualified Synthesizer.Plain.Filter.Recursive.Allpass as Allpass +-- import qualified Synthesizer.Plain.Signal as Sig + +-- import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty + +import Test.QuickCheck (quickCheck, {- Property, (==>) -}) +import Test.Utility (equalList, ) + +-- import qualified Algebra.Module as Module +-- import qualified Algebra.RealField as RealField +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive + +import Control.Monad.Trans.State (runState, ) + +-- import Debug.Trace (trace, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +{- this will not work due to the poles +parameter :: Double -> Double -> Bool +parameter phase freq = + approxEqual eps phase + (Allpass.makePhase (Allpass.parameter phase freq) freq) +-} + + +cascadeStep :: Rational -> Rational -> (Rational, Rational, [Rational]) -> Bool +cascadeStep k u (s0,s1,ns) = + let p = Allpass.Parameter k + s = s0:s1:ns + in equalList $ + runState (Allpass.cascadeStepStack p u) s : + runState (Allpass.cascadeStepRec p u) s : + runState (Allpass.cascadeStepScanl p u) s : + [] + + +cascade :: NonNeg.Int -> Sig.T Rational -> Sig.T Rational -> Bool +cascade order ks xs = + let ps = map Allpass.Parameter ks + n = NonNeg.toNumber order + in Allpass.cascadeState n ps xs == + Allpass.cascadeIterative n ps xs + + +tests :: [(String, IO ())] +tests = + ("cascadeStep", quickCheck cascadeStep) : + ("cascade", quickCheck cascade) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs b/test/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Filter/Hilbert.hs @@ -0,0 +1,44 @@ +module Test.Sound.Synthesizer.Plain.Filter.Hilbert (tests) where + +import qualified Synthesizer.Plain.Filter.Recursive.Hilbert as Hilbert +import qualified Synthesizer.Plain.Filter.Recursive.Allpass as Allpass +import qualified Synthesizer.Plain.Signal as Sig + +import qualified Synthesizer.Causal.Process as Causal + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty + +import Test.QuickCheck (quickCheck, {- Property, (==>) -}) +-- import Test.Utility (equalList, ) + +-- import qualified Algebra.Module as Module +-- import qualified Algebra.RealField as RealField +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive +-- import qualified Number.Complex as Complex + +import Data.Tuple.HT (mapPair, ) + +-- import Debug.Trace (trace, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +cascade :: NonEmpty.T (Rational, Rational) -> Sig.T Rational -> Bool +cascade ks xs = + let p = uncurry Hilbert.Parameter $ unzip $ + map (mapPair (Allpass.Parameter, Allpass.Parameter)) $ + NonEmpty.toList ks + in Hilbert.run2 p xs == + Causal.apply (Hilbert.causal2 p) xs +{- + in map Complex.real (Hilbert.run2 p xs) == xs +-} + + +tests :: [(String, IO ())] +tests = + ("hilbert", quickCheck cascade) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Interpolation.hs b/test/Test/Sound/Synthesizer/Plain/Interpolation.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Interpolation.hs @@ -0,0 +1,343 @@ +module Test.Sound.Synthesizer.Plain.Interpolation ( + T, ip, + LinePreserving, lpIp, + tests, + use, useLP, use2, + -- only for debugging + frequencyModulationBackCompare, + frequencyModulationForth0Compare, + frequencyModulationStorableChunkSizeCompare, + frequencyModulationStorableCompare, + ) where + +import qualified Synthesizer.Plain.Interpolation as Interpolation +import qualified Synthesizer.Interpolation.Class as Interpol +import qualified Synthesizer.Interpolation.Custom as ExampleCustom +import qualified Synthesizer.Interpolation.Module as ExampleModule +import qualified Synthesizer.Interpolation as InterpolationCore + +import qualified Synthesizer.Causal.Interpolation as InterpolC +import qualified Synthesizer.Causal.Process as Causal +import qualified Synthesizer.Generic.Filter.NonRecursive as FiltG +import qualified Synthesizer.Generic.Signal as SigG +import qualified Synthesizer.State.Filter.NonRecursive as FiltS +import qualified Synthesizer.State.Signal as SigS + +import qualified Synthesizer.Storable.Filter.NonRecursive as FiltSt +import qualified Synthesizer.Storable.Signal as SigSt + +import Test.QuickCheck (quickCheck, Arbitrary(arbitrary), elements, {- Property, (==>), -} Testable, ) +-- import Test.Utility + +import Foreign.Storable (Storable, ) + +import qualified Algebra.VectorSpace as VectorSpace +import qualified Algebra.Module as Module +import qualified Algebra.RealField as RealField +import qualified Algebra.Field as Field +import qualified Algebra.RealRing as RealRing +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty +import qualified Data.List.Match as Match +import Control.Monad (liftM2, ) + +import Test.Utility (equalList, ArbChar, unpackArbString, ) + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + + +instance Arbitrary InterpolationCore.Margin where + arbitrary = + liftM2 InterpolationCore.Margin + (fmap abs arbitrary) + (fmap abs arbitrary) + + +use :: + (Interpolation.T a v -> x) -> + (T a v -> x) +use f ipt = + f (ip ipt) + +useLP :: + (Interpolation.T a v -> x) -> + (LinePreserving a v -> x) +useLP f ipt = + f (lpIp ipt) + +use2 :: + (Interpolation.T a v -> + Interpolation.T a v -> x) -> + (T a v -> + T a v -> x) +use2 f = + use $ \ ipLeap -> + use $ \ ipStep -> + f ipLeap ipStep + + + +data T a v = Cons {name :: String, ip :: Interpolation.T a v} + +instance Show (T a v) where + show x = name x + +instance (Field.C a, Interpol.C a v) => Arbitrary (T a v) where + arbitrary = elements $ + Cons "constant" ExampleCustom.constant : + Cons "linear" ExampleCustom.linear : + Cons "cubic" ExampleCustom.cubic : + [] + + + +data LinePreserving a v = + LPCons {lpName :: String, lpIp :: Interpolation.T a v} + +instance Show (LinePreserving a v) where + show x = lpName x + +instance (Field.C a, Interpol.C a v) => Arbitrary (LinePreserving a v) where + arbitrary = elements $ + LPCons "linear" ExampleCustom.linear : + LPCons "cubic" ExampleCustom.cubic : + [] + + + +constant :: + (Interpol.C a v, Module.C a v, Eq v) => + a -> v -> [v] -> Bool +constant t x0 xs = + equalList $ map ($(x0:xs)) $ map ($t) $ + Interpolation.func ExampleCustom.constant : + Interpolation.func ExampleCustom.piecewiseConstant : + Interpolation.func ExampleModule.constant : + Interpolation.func ExampleModule.piecewiseConstant : + [] + +linear :: + (Interpol.C a v, Module.C a v, Eq v) => + a -> v -> v -> [v] -> Bool +linear t x0 x1 xs = + equalList $ map ($(x0:x1:xs)) $ map ($t) $ + Interpolation.func ExampleCustom.linear : + Interpolation.func ExampleCustom.piecewiseLinear : + Interpolation.func ExampleModule.linear : + Interpolation.func ExampleModule.piecewiseLinear : + [] + +cubic :: + (Interpol.C a v, VectorSpace.C a v, Eq v) => + a -> v -> v -> v -> v -> [v] -> Bool +cubic t x0 x1 x2 x3 xs = + equalList $ map ($(x0:x1:x2:x3:xs)) $ map ($t) $ + Interpolation.func ExampleCustom.cubic : + Interpolation.func ExampleCustom.piecewiseCubic : + Interpolation.func ExampleModule.cubic : + Interpolation.func ExampleModule.cubicAlt : + Interpolation.func ExampleModule.piecewiseCubic : + [] + + +controlAboveOne :: (RealRing.C t) => [t] -> [t] +controlAboveOne = + map ((one+) . abs) + +frequencyModulationForth0 :: + (RealField.C t, Eq v) => + [t] -> [v] -> Bool +frequencyModulationForth0 cs0 xs = + let cs = controlAboveOne cs0 + in Causal.apply + (InterpolC.relative ExampleModule.constant zero + (FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) (SigS.fromList xs))) + (Match.take xs cs) + == Match.take cs xs + +frequencyModulationForth0Compare :: + (RealField.C t, Eq v) => + [t] -> [v] -> ([v], [v], [v]) +frequencyModulationForth0Compare cs0 xs = + let cs = controlAboveOne cs0 + in (Match.take cs + (Causal.apply + (InterpolC.relative ExampleModule.constant zero + (FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) (SigS.fromList xs))) + (Match.take xs cs)), + SigS.toList + (FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) (SigS.fromList xs)), + Match.take cs xs) + + +frequencyModulationForth1 :: + (RealField.C t, Eq v) => + [t] -> [v] -> Bool +frequencyModulationForth1 cs0 xs = + case controlAboveOne cs0 of + [] -> True + (c:cs) -> + Causal.apply + (InterpolC.relative ExampleModule.constant c + (FiltS.inverseFrequencyModulationFloor + (SigS.fromList ((c+one):cs)) (SigS.fromList xs))) + (Match.take xs cs) + == Match.take cs xs + + + +controlBelowOne :: (RealField.C t) => [t] -> [t] +controlBelowOne = + map fraction + + +frequencyModulationBack :: + (RealField.C t, Eq v) => + [t] -> NonEmpty.T v -> Bool +frequencyModulationBack cs0 xs0 = + let cs = controlBelowOne cs0 + xs = NonEmpty.toInfiniteList xs0 + in take (floor (sum cs)) xs == + (SigS.toList $ + FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) + (SigS.fromList $ + Causal.apply + (InterpolC.relative ExampleModule.constant zero + (SigS.fromList xs)) + cs)) + + +frequencyModulationBackCompare :: + (RealField.C t, Eq v) => + [t] -> [v] -> (SigS.T v, SigS.T v) +frequencyModulationBackCompare cs0 xs = + let cs = controlBelowOne cs0 + in (FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) + (SigS.fromList $ + Causal.apply + (InterpolC.relative ExampleModule.constant zero + (SigS.fromList (cycle xs))) + cs), + SigS.fromList $ + Causal.apply + (InterpolC.relative ExampleModule.constant zero + (SigS.fromList (cycle xs))) + cs) + +frequencyModulationGeneric :: + (RealField.C t, Eq v) => + [t] -> [v] -> Bool +frequencyModulationGeneric cs xs = + SigS.toList + (FiltS.inverseFrequencyModulationFloor + (SigS.fromList cs) (SigS.fromList xs)) + == FiltG.inverseFrequencyModulationFloor + SigG.defaultLazySize cs xs + + +makeChunkSize :: Int -> SigSt.ChunkSize +makeChunkSize size = + SigSt.chunkSize (1 + abs size) + +{- +makeExactFraction :: (Int,Int) -> Double +makeExactFraction (n,d) = + fromIntegral n * 2 ^- (- mod (fromIntegral d) 4) +-} + +frequencyModulationStorableChunkSize :: + (Storable v, RealField.C t, Eq v) => + Int -> Int -> + Int -> Int -> + [t] -> [v] -> + Bool +frequencyModulationStorableChunkSize size0 size1 xsize0 xsize1 cs xs = + FiltSt.inverseFrequencyModulationFloor + (makeChunkSize size0) cs + (SigSt.fromList (makeChunkSize xsize0) xs) + == + FiltSt.inverseFrequencyModulationFloor + (makeChunkSize size1) cs + (SigSt.fromList (makeChunkSize xsize1) xs) + + +frequencyModulationStorableChunkSizeCompare :: + (Storable v, RealField.C t, Eq v) => + Int -> Int -> + Int -> Int -> + [t] -> [v] -> + (SigSt.T v, SigSt.T v) +frequencyModulationStorableChunkSizeCompare size0 size1 xsize0 xsize1 cs xs = + (FiltSt.inverseFrequencyModulationFloor + (makeChunkSize size0) cs + (SigSt.fromList (makeChunkSize xsize0) xs), + FiltSt.inverseFrequencyModulationFloor + (makeChunkSize size1) cs + (SigSt.fromList (makeChunkSize xsize1) xs)) + + +frequencyModulationStorable :: + (Storable v, RealField.C t, Eq v) => + Int -> Int -> + [t] -> [v] -> + Bool +frequencyModulationStorable size xsize cs xs = + SigSt.toList + (FiltSt.inverseFrequencyModulationFloor (makeChunkSize size) cs + (SigSt.fromList (makeChunkSize xsize) xs)) + == FiltG.inverseFrequencyModulationFloor + SigG.defaultLazySize cs xs + + +frequencyModulationStorableCompare :: + (Storable v, RealField.C t, Eq v) => + Int -> Int -> + [t] -> [v] -> + ([v], SigSt.T v) +frequencyModulationStorableCompare size xsize cs xs = + (FiltG.inverseFrequencyModulationFloor + SigG.defaultLazySize cs xs, + FiltSt.inverseFrequencyModulationFloor (makeChunkSize size) cs + (SigSt.fromList (makeChunkSize xsize) xs)) + + + +testRational :: + (Testable t) => + (Rational -> Rational -> t) -> IO () +testRational = quickCheck + +testFM :: + (Testable t, Arbitrary (sigX ArbChar), Show (sigX ArbChar)) => + ([Rational] -> sigX ArbChar -> t) -> IO () +testFM = quickCheck + +tests :: [(String, IO ())] +tests = + ("constant", testRational constant) : + ("linear", testRational linear ) : + ("cubic", testRational cubic ) : + ("frequencyModulationForth0", testFM frequencyModulationForth0) : + ("frequencyModulationForth1", testFM frequencyModulationForth1) : + ("frequencyModulationBack", testFM frequencyModulationBack) : + ("frequencyModulationGeneric", testFM frequencyModulationGeneric) : + ("frequencyModulationStorableChunkSize", + quickCheck (\size0 size1 xsize0 xsize1 cs xs -> + frequencyModulationStorableChunkSize size0 size1 xsize0 xsize1 + (cs::[Rational]) (unpackArbString xs))) : + ("frequencyModulationStorable", + quickCheck (\size xsize cs xs -> + frequencyModulationStorable size xsize + (cs::[Rational]) (unpackArbString xs))) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/NonEmpty.hs b/test/Test/Sound/Synthesizer/Plain/NonEmpty.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/NonEmpty.hs @@ -0,0 +1,34 @@ +module Test.Sound.Synthesizer.Plain.NonEmpty where + +import Test.QuickCheck (Arbitrary, arbitrary, ) +import Control.Monad (liftM2, ) + + +data T a = Cons a [a] + +toList :: T a -> [a] +toList (Cons x xs) = + (x:xs) + +toInfiniteList :: T a -> [a] +toInfiniteList = + cycle . toList + +instance Functor T where + fmap f (Cons x xs) = + Cons (f x) (map f xs) + +instance Arbitrary a => Arbitrary (T a) where + arbitrary = liftM2 Cons arbitrary arbitrary + +instance Show a => Show (T a) where + showsPrec p (Cons x xs) = + showsPrec p (x:xs) + +{- +instance Show a => Show (T a) where + showsPrec p (Cons x xs) = + showParen (p >= 10) $ + showString "cycle " . + showsPrec 11 (x:xs) +-} diff --git a/test/Test/Sound/Synthesizer/Plain/Oscillator.hs b/test/Test/Sound/Synthesizer/Plain/Oscillator.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Oscillator.hs @@ -0,0 +1,39 @@ +module Test.Sound.Synthesizer.Plain.Oscillator (tests) where + +import qualified Synthesizer.Plain.Oscillator as Osci +import qualified Synthesizer.Basic.Wave as Wave +-- import qualified Synthesizer.Plain.Interpolation as Interpolation + +import qualified Test.Sound.Synthesizer.Plain.Wave as WaveTest +-- import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest + +import Test.QuickCheck (quickCheck, {- Property, (==>), -} ) + +import qualified Algebra.RealField as RealField + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + + +phaseShapeMod :: (RealField.C a, Eq b) => (Wave.T a b) -> a -> [a] -> Bool +phaseShapeMod wave freq phases = + Osci.phaseMod wave freq phases == + Osci.shapeMod (Wave.phaseOffset wave) zero freq phases + +phaseShapeModRational :: + WaveTest.Ring Rational -> Integer -> Integer -> [Integer] -> Bool +phaseShapeModRational w denom0 freq0 phases0 = + let denom = 1 + abs denom0 + freq = freq0 % denom + phases = map (% denom) phases0 + in phaseShapeMod (WaveTest.ringWave w) freq phases + + + +tests :: [(String, IO ())] +tests = + ("phaseShapeModRational", quickCheck phaseShapeModRational) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/ToneModulation.hs b/test/Test/Sound/Synthesizer/Plain/ToneModulation.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/ToneModulation.hs @@ -0,0 +1,478 @@ +module Test.Sound.Synthesizer.Plain.ToneModulation (tests, ) where + +import Test.Sound.Synthesizer.Basic.ToneModulation ( + minLength, + minLengthMargin, + shapeLimits, + testRationalLineIp, + testRationalIp, + ) + +import qualified Synthesizer.Plain.Oscillator as Osci +import qualified Synthesizer.Plain.Interpolation as Interpolation +import qualified Synthesizer.Plain.ToneModulation as ToneModL +import qualified Synthesizer.Plain.Wave as WaveL +import Synthesizer.Interpolation (marginNumber, ) + +import qualified Synthesizer.Basic.Wave as Wave +import qualified Synthesizer.Basic.Phase as Phase + +import qualified Test.Sound.Synthesizer.Plain.NonEmpty as NonEmpty +import qualified Test.Sound.Synthesizer.Plain.Interpolation as InterpolationTest + +import Test.QuickCheck (quickCheck, Property, (==>), ) +import Test.Utility (ArbChar, ) + +import qualified Number.NonNegative as NonNeg +import qualified Number.NonNegativeChunky as Chunky + +import qualified Algebra.RealTranscendental as RealTrans +import qualified Algebra.Module as Module +import qualified Algebra.RealField as RealField +import qualified Algebra.Additive as Additive +import qualified Algebra.ZeroTestable as ZeroTestable + +import Data.List.HT (isAscending, ) +import Data.Ord.HT (limit, ) +import Data.Tuple.HT (mapPair, mapSnd, ) +import qualified Data.List as List + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +{- +Properties that do not hold: + commutativity of limitRelativeShapes and integrateFractional: + Does not hold because when you clip the integral skips at the end, + you would have to clear the fractional part, too. +-} + + + +absolutize :: (Additive.C a) => a -> [a] -> [a] +absolutize = scanl (+) + +limitMinRelativeValues :: + Int -> Int -> [NonNeg.Int] -> Bool +limitMinRelativeValues xMin x0 xsnn = + let xs = map NonNeg.toNumber xsnn + in map (max xMin) (absolutize x0 xs) == + uncurry absolutize (ToneModL.limitMinRelativeValues xMin x0 xs) + +limitMaxRelativeValues :: + Int -> Int -> [NonNeg.Int] -> Bool +limitMaxRelativeValues xMax x0 xsnn = + let xs = map NonNeg.toNumber xsnn + in map (min xMax) (absolutize x0 xs) == + uncurry absolutize (ToneModL.limitMaxRelativeValues xMax x0 xs) + +limitMaxRelativeValuesNonNeg :: + Int -> Int -> [NonNeg.Int] -> Bool +limitMaxRelativeValuesNonNeg xMax x0 xsnn = + let xs = map NonNeg.toNumber xsnn + in map (min xMax) (absolutize x0 xs) == + uncurry absolutize (ToneModL.limitMaxRelativeValuesNonNeg xMax x0 xs) + +-- chunky type is not necessary here but testing it a little is not wrong +limitMinRelativeValuesIdentity :: + Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool +limitMinRelativeValuesIdentity x0 xs = + (x0,xs) == ToneModL.limitMinRelativeValues 0 x0 xs + +limitMaxRelativeValuesIdentity :: + Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool +limitMaxRelativeValuesIdentity x0 xs = + let inf = 1 + inf + in (x0,xs) == ToneModL.limitMaxRelativeValues inf x0 xs + +limitMaxRelativeValuesNonNegIdentity :: + Chunky.T NonNeg.Int -> [Chunky.T NonNeg.Int] -> Bool +limitMaxRelativeValuesNonNegIdentity x0 xs = + let inf = 1 + inf + in (x0,xs) == ToneModL.limitMaxRelativeValuesNonNeg inf x0 xs + +limitMaxRelativeValuesInfinity :: + Chunky.T NonNeg.Int -> NonEmpty.T (Chunky.T NonNeg.Int) -> Bool +limitMaxRelativeValuesInfinity x0 ixs = + let inf = 1 + inf + ys = NonEmpty.toInfiniteList ixs + (z0,zs) = ToneModL.limitMaxRelativeValues inf x0 ys + in (x0, take 100 ys) == (z0, take 100 zs) + +limitMaxRelativeValuesNonNegInfinity :: + Chunky.T NonNeg.Int -> NonEmpty.T (Chunky.T NonNeg.Int) -> Bool +limitMaxRelativeValuesNonNegInfinity x0 ixs = + let inf = 1 + inf + ys = NonEmpty.toInfiniteList ixs + (z0,zs) = ToneModL.limitMaxRelativeValuesNonNeg inf x0 ys + in (x0, take 100 ys) == (z0, take 100 zs) + + +dropRem :: Eq a => NonNeg.Int -> [a] -> Bool +dropRem nn xs = + let n = NonNeg.toNumber nn + in map (flip ToneModL.dropRem xs) [0 .. n + length xs] == + map ((,) 0) (List.tails xs) ++ map (flip (,) []) [1..n] + + +sampledToneSine :: (RealTrans.C a, Module.C a a) => + NonNeg.T a -> NonNeg.Int -> a -> a -> a -> Bool +sampledToneSine periodNN ext phase0 shape phase = + let ipLeap = Interpolation.cubic + ipStep = Interpolation.cubic + ten = fromInteger 10 + period = ten + NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (Osci.staticSine phase0 (recip period)) + in abs (WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` (Phase.fromRepresentative phase) - + head (Osci.staticSine (phase0+phase) zero)) < ten ^- (-2) + + +sampledToneSineList :: (RealTrans.C a, Module.C a a) => + NonNeg.T a -> NonNeg.Int -> a -> a -> [a] -> [a] -> Bool +sampledToneSineList periodNN ext origPhase phase shapes freqs = + let ipLeap = Interpolation.cubic + ipStep = Interpolation.cubic + ten = fromInteger 10 + period = ten + NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (Osci.staticSine origPhase (recip period)) + in all ((< ten ^- (-2)) . abs) $ + zipWith (-) + (Osci.shapeFreqMod (WaveL.sampledTone ipLeap ipStep period tone) + phase shapes freqs) + (Osci.freqModSine (origPhase+phase) freqs) + + +sampledToneLinear :: (RealField.C a, Module.C a v, Eq v) => + InterpolationTest.LinePreserving a v -> + InterpolationTest.LinePreserving a v -> + NonNeg.T a -> NonNeg.Int -> (v,v) -> a -> Phase.T a -> Property +sampledToneLinear = + InterpolationTest.useLP $ \ ipLeap -> + InterpolationTest.useLP $ \ ipStep -> + \ periodNN ext (i,d) shape phase -> + let period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + ramp = take len (List.iterate (d+) i) + limits = + mapPair (fromIntegral, fromIntegral) $ + shapeLimits ipLeap ipStep periodInt len + in period /= zero ==> + -- should be (fraction phase), right? + WaveL.sampledTone ipLeap ipStep period ramp shape `Wave.apply` phase == + i + limit limits shape *> d +{- +let len=100; period=1/0.06::Double; ip = Interpolation.linear in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (0,fromIntegral len)) [\s -> WaveL.sampledTone ip ip period (take len $ iterate (1+) (0::Double)) s 0, limit (mapPair (fromIntegral, fromIntegral) $ shapeLimits ip ip (round period::Int) len)] +-} + +sampledToneStair :: (RealField.C a, Module.C a v, Eq v) => + InterpolationTest.LinePreserving a v -> + NonNeg.Int -> NonNeg.Int -> (v,v) -> a -> Property +sampledToneStair = + InterpolationTest.useLP $ \ ipLeap + periodIntNN ext (i,d) shape -> + let ipStep = Interpolation.constant + periodInt = NonNeg.toNumber periodIntNN + period = fromIntegral periodInt + len0 = minLength ipLeap ipStep periodInt ext + (rep,rm) = divMod (negate len0) periodInt + len = len0 + rm + stair = + concatMap (replicate periodInt) $ + take (negate rep) (List.iterate (period*>d+) i) + limits = + mapPair (fromIntegral, fromIntegral) $ + shapeLimits ipLeap ipStep periodInt len + in periodInt /= zero ==> + WaveL.sampledTone ipLeap ipStep period stair shape `Wave.apply` zero == + i + limit limits shape *> d +{- +let len=periodInt*rep; rep=10; periodInt = 14::Int; period=fromIntegral periodInt; ipl = Interpolation.linear; ipc = Interpolation.constant in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (-10,10+fromIntegral len)) [\s -> WaveL.sampledTone ipl ipc period (concatMap (replicate periodInt) $ take rep $ iterate (period+) (0::Double)) s 0, limit (mapPair (fromIntegral, fromIntegral) $ shapeLimits ipl ipc periodInt len)] +-} + +{- +sampledToneSaw :: (RealField.C a, Module.C a v, Eq v) => + InterpolationTest.LinePreserving a v -> + InterpolationTest.T a v -> + NonNeg.Int -> NonNeg.Int -> (v,v) -> a -> a -> Property +sampledToneSaw iptLeap iptStep periodIntNN ext (i,d) shape phase = + let ipLeap = InterpolationTest.lpIp iptLeap + ipStep = InterpolationTest.ip iptStep + periodInt = NonNeg.toNumber periodIntNN + period = fromIntegral periodInt + len0 = minLength ipLeap ipStep periodInt ext + rep = negate $ div (negate len0) periodInt + saw = + concat $ replicate rep $ + take periodInt $ List.iterate (d+) i + in periodInt /= zero ==> + WaveL.sampledTone ipLeap ipStep period saw shape phase == + i + fraction phase *> d +-} + +sampledToneStatic :: (RealField.C a, Eq v) => + InterpolationTest.T a v -> + InterpolationTest.T a v -> + NonNeg.Int -> (v,[v]) -> a -> a -> Property +sampledToneStatic = + InterpolationTest.use2 $ \ ipLeap ipStep + ext (x,xs) shape phase -> + let wave = x:xs + periodInt = length wave + period = fromIntegral periodInt + len = minLength ipLeap ipStep periodInt ext + rep = negate $ div (negate len) periodInt + tone = concat $ replicate rep wave + in period /= zero ==> + WaveL.sampledTone ipLeap ipStep period tone shape `Wave.apply` (Phase.fromRepresentative phase) == + Interpolation.cyclicPad Interpolation.single ipStep (phase*period) wave +{- +let wave = [1,-1,0.5,-0.5::Double]; period = fromIntegral (length wave) :: Double; ip = Interpolation.linear in GNUPlot.plotFuncs [] (GNUPlot.linearScale 1000 (-1,3)) [WaveL.sampledTone ip ip period (concat $ replicate 3 wave) 0.3, \phase -> Interpolation.cyclicPad Interpolation.single Interpolation.linear (phase*period) wave] +-} + + + +shapeFreqModFromSampledToneLimitIdentity :: (RealField.C t) => + Interpolation.Margin -> + Interpolation.Margin -> + NonNeg.Int -> NonEmpty.T y -> (t, NonEmpty.T (NonNeg.T t)) -> Bool +shapeFreqModFromSampledToneLimitIdentity + marginLeap marginStep periodIntNN ixs (shape0,shapesNN) = + let periodInt = NonNeg.toNumber periodIntNN + shapes = fmap NonNeg.toNumber shapesNN + a = snd + (ToneModL.limitRelativeShapes + marginLeap marginStep + periodInt (NonEmpty.toInfiniteList ixs) + (shape0, NonEmpty.toInfiniteList shapes)) !! 100 + in a == a + + +oscillatorCoords :: (RealField.C t) => + NonNeg.Int -> NonNeg.T t -> t -> Phase.T t -> [NonNeg.T t] -> [t] -> Property +oscillatorCoords + periodIntNN periodNN shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + periodRound = fromIntegral periodInt + coords = + ToneModL.oscillatorCoords + periodInt period + (shape0, shapes) (phase, freqs) + in period /= zero && periodInt /= zero ==> + all + (\(skip,(k,(qShape,qWave))) -> + skip >= zero && + isAscending [negate periodInt, k, zero] && + isAscending [zero, qShape, one] && + isAscending [zero, qWave, periodRound]) + (tail coords) + + +shapeFreqModFromSampledToneCoordsIdentity :: + (RealField.C t, ZeroTestable.C t) => + NonNeg.Int -> NonNeg.T t -> (t, [NonNeg.T t]) -> Property +shapeFreqModFromSampledToneCoordsIdentity + periodIntNN periodNN (shape0,shapesNN) = + let period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + shapes = map NonNeg.toNumber shapesNN + phase = Phase.fromRepresentative $ shape0 / period + freqs = map (/period) shapes + in period /= zero ==> + all + (isZero . fst . snd . snd) + (ToneModL.oscillatorCoords + periodInt period (shape0, shapes) (phase, freqs)) + + +shapeFreqModFromSampledTone :: (RealField.C t, Eq v) => + InterpolationTest.T t v -> + InterpolationTest.T t v -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> + Property +shapeFreqModFromSampledTone = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs shape0 phase shapesNN freqs -> + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledToneA = + Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone + shape0 phase shapes freqs + resampledToneB = + Osci.shapeFreqMod + (WaveL.sampledTone ipLeap ipStep period tone) + phase (scanl (+) shape0 shapes) freqs + in period /= zero ==> + resampledToneA == resampledToneB +{- +let len=100; period=1/0.06::Double; ip = Interpolation.linear; tone = take len $ iterate (1+) (0::Double); shape0=0; shapes = replicate 100 1; in GNUPlot.plotLists [] [Osci.shapeFreqMod (WaveL.sampledTone ip ip period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ip ip period tone shape0 0 shapes (repeat 0)] +*Test.Sound.Synthesizer.Plain.Oscillator> let len=100; period=1/0.06::Double; ip = Interpolation.linear; tone = take len $ iterate (1+) (0::Double); shape0=0; shapes = concat $ replicate 50 [1.5,0.5]; in GNUPlot.plotLists [] [Osci.shapeFreqMod (WaveL.sampledTone ip ip period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ip ip period tone shape0 0 shapes (repeat 0)] +*Test.Sound.Synthesizer.Plain.Oscillator> let len=100; period=1/0.06::Rational; ipLeap = Interpolation.linear; ipStep = Interpolation.constant; tone = take len $ iterate (1+) (0::Rational); shape0=0; shapes = concat $ replicate 50 [1.5,0.5]; in GNUPlot.plotLists [] (map (map (\x -> fromRational' x :: Double)) [Osci.shapeFreqMod (WaveL.sampledTone ipLeap ipStep period tone) 0 (scanl (+) shape0 shapes) (repeat 0), Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone shape0 0 shapes (repeat 0)]) +-} + + +shapePhaseFreqModFromSampledTone :: (RealField.C t, Eq v) => + InterpolationTest.T t v -> + InterpolationTest.T t v -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> [t] -> + Property +shapePhaseFreqModFromSampledTone = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs shape0 phase shapesNN phaseDistorts freqs -> + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + resampledToneA = + Osci.shapePhaseFreqModFromSampledTone ipLeap ipStep period tone + shape0 phase shapes phaseDistorts freqs + resampledToneB = + Osci.shapeFreqMod + (uncurry $ + Wave.phaseOffset . + WaveL.sampledTone ipLeap ipStep period tone) + phase (zip (scanl (+) shape0 shapes) phaseDistorts) freqs + in period /= zero ==> + resampledToneA == resampledToneB + + +oscillatorCells :: (RealField.C t, Eq v) => + Interpolation.Margin -> + Interpolation.Margin -> + NonNeg.Int -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + t -> t -> [NonNeg.T t] -> [t] -> + Property +oscillatorCells + marginLeap marginStep periodIntNN periodNN ext ixs shape0 phase shapesNN freqs = + let shapes = map NonNeg.toNumber shapesNN + period = NonNeg.toNumber periodNN + periodInt = NonNeg.toNumber periodIntNN + len = minLengthMargin marginLeap marginStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + crop = cropCell marginLeap marginStep + resampledToneA = + ToneModL.oscillatorCells + marginLeap marginStep periodInt period tone + (shape0, shapes) (Phase.fromRepresentative phase, freqs) + resampledToneB = + Osci.shapeFreqMod + (Wave.Cons . ToneModL.sampledToneCell + (ToneModL.makePrototype marginLeap marginStep + periodInt period tone)) + phase (scanl (+) shape0 shapes) freqs + in period /= zero && + periodInt /= zero && + marginNumber marginLeap > zero && + marginNumber marginStep > zero ==> + map crop resampledToneA == map crop resampledToneB + +cropCell :: + Interpolation.Margin -> + Interpolation.Margin -> + ((t,t), ToneModL.Cell v) -> ((t,t), ToneModL.Cell v) +cropCell ipLeap ipStep = + mapSnd + (take (marginNumber ipStep) . + map (take (marginNumber ipLeap))) + + +shapeFreqModFromSampledToneIdentity :: (RealField.C t, Eq v) => + InterpolationTest.T t v -> + InterpolationTest.T t v -> + NonNeg.T t -> + NonNeg.Int -> NonEmpty.T v -> + Property +shapeFreqModFromSampledToneIdentity = + InterpolationTest.use2 $ \ ipLeap ipStep + periodNN ext ixs -> + let period = NonNeg.toNumber periodNN + periodInt = round period + len = minLength ipLeap ipStep periodInt ext + tone = take len (NonEmpty.toInfiniteList ixs) + shape0 = zero + shapes = repeat one + phase = zero + freqs = repeat (recip period) + (n0,n1) = + shapeLimits ipLeap ipStep periodInt len + + resampledTone = + Osci.shapeFreqModFromSampledTone ipLeap ipStep period tone + shape0 phase shapes freqs + in period /= zero ==> + and (drop n0 (take (succ n1) (zipWith (==) resampledTone tone))) + + +tests :: [(String, IO ())] +tests = + ("limitMinRelativeValues", quickCheck limitMinRelativeValues) : + ("limitMaxRelativeValues", quickCheck limitMaxRelativeValues) : + ("limitMaxRelativeValuesNonNeg", + quickCheck limitMaxRelativeValuesNonNeg) : + ("limitMinRelativeValuesIdentity", + quickCheck limitMinRelativeValuesIdentity) : + ("limitMaxRelativeValuesIdentity", + quickCheck limitMaxRelativeValuesIdentity) : + ("limitMaxRelativeValuesNonNegIdentity", + quickCheck limitMaxRelativeValuesNonNegIdentity) : + ("limitMaxRelativeValuesInfinity", + quickCheck limitMaxRelativeValuesInfinity) : + ("limitMaxRelativeValuesNonNegInfinity", + quickCheck limitMaxRelativeValuesNonNegInfinity) : + ("dropRem", quickCheck (dropRem :: NonNeg.Int -> [ArbChar] -> Bool)) : + ("sampledToneSine", + quickCheck (\period -> sampledToneSine (period :: NonNeg.Double))) : + ("sampledToneSineList", + quickCheck (\period -> sampledToneSineList (period :: NonNeg.Double))) : + ("sampledToneLinear", + testRationalLineIp sampledToneLinear) : + ("sampledToneStair", + testRationalLineIp sampledToneStair) : +{- + ("sampledToneSaw", + testRationalLineIp sampledToneSaw) : +-} + ("sampledToneStatic", + testRationalIp sampledToneStatic) : + ("shapeFreqModFromSampledToneLimitIdentity", + quickCheck (\ml ms p ixs (t,ts) -> + shapeFreqModFromSampledToneLimitIdentity ml ms p + (ixs::NonEmpty.T Rational) (t::Rational,ts))) : + ("oscillatorCoords", + quickCheck (\periodInt period -> + oscillatorCoords + periodInt (period :: NonNeg.Rational))) : + ("shapeFreqModFromSampledToneCoordsIdentity", + quickCheck (\periodInt period -> + shapeFreqModFromSampledToneCoordsIdentity + periodInt (period :: NonNeg.Rational))) : + ("shapeFreqModFromSampledTone", + testRationalIp shapeFreqModFromSampledTone) : + ("shapePhaseFreqModFromSampledTone", + testRationalIp shapePhaseFreqModFromSampledTone) : + ("oscillatorCells", + quickCheck (\ml ms periodInt period ext ixs -> + oscillatorCells ml ms periodInt (period :: NonNeg.Rational) + ext (ixs :: NonEmpty.T ArbChar))) : + ("shapeFreqModFromSampledToneIdentity", + testRationalIp shapeFreqModFromSampledToneIdentity) : + [] diff --git a/test/Test/Sound/Synthesizer/Plain/Wave.hs b/test/Test/Sound/Synthesizer/Plain/Wave.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Plain/Wave.hs @@ -0,0 +1,75 @@ +module Test.Sound.Synthesizer.Plain.Wave (Ring, ringWave, tests) where + +import qualified Synthesizer.Basic.Wave as Wave +import qualified Synthesizer.Basic.Phase as Phase + +import Test.QuickCheck (quickCheck, Arbitrary(arbitrary), elements, oneof, choose, {- Property, (==>), -} ) +-- import Test.Utility + +import qualified Number.NonNegative as NonNeg + +import qualified Algebra.RealTranscendental as RealTrans +import qualified Algebra.Ring as Ring + +import Control.Monad (liftM, liftM2, ) +import System.Random (Random) + + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + + + +data Ring a = Ring {ringName :: String, ringWave :: Wave.T a a} + +instance Show (Ring a) where + show = ringName + +instance (Ord a, Ring.C a) => Arbitrary (Ring a) where + arbitrary = elements $ + Ring "saw" Wave.saw : + Ring "square" Wave.square : + Ring "triangle" Wave.triangle : + [] + + + + +data ZeroDCOffset a = ZeroDCOffset {zdcName :: String, zdcWave :: Wave.T a a} + +instance Show (ZeroDCOffset a) where + show = zdcName + +instance (RealTrans.C a, Random a) => Arbitrary (ZeroDCOffset a) where + arbitrary = + let cons n w = return (ZeroDCOffset n w) + in oneof $ + cons "sine" Wave.sine : + cons "saw" Wave.saw : + cons "square" Wave.square : + cons "triangle" Wave.triangle : + liftM + (ZeroDCOffset "squareBalanced" . Wave.squareBalanced) + (choose (negate one, one)) : + liftM2 + (\w r -> ZeroDCOffset "trapezoidBalanced" (Wave.trapezoidBalanced w r)) + (choose (zero, one)) + (choose (negate one, one)) : + [] + + +zeroDCOffset :: ZeroDCOffset Double -> NonNeg.Int -> Bool +zeroDCOffset w periodIntNN = + let periodInt = 100 + NonNeg.toNumber periodIntNN + period = fromIntegral periodInt + xs = take periodInt $ map Phase.fromRepresentative $ + map (/period) $ iterate (1+) 0.5 + in abs (sum (map (Wave.apply (zdcWave w)) xs)) < period / fromInteger 100 + + +tests :: [(String, IO ())] +tests = + ("zeroDCOffset", quickCheck zeroDCOffset) : + [] diff --git a/test/Test/Sound/Synthesizer/Storable/Cut.hs b/test/Test/Sound/Synthesizer/Storable/Cut.hs new file mode 100644 --- /dev/null +++ b/test/Test/Sound/Synthesizer/Storable/Cut.hs @@ -0,0 +1,40 @@ +module Test.Sound.Synthesizer.Storable.Cut (tests) where + +import qualified Synthesizer.Storable.Cut as CutSt +import qualified Synthesizer.Storable.Signal as SigSt + +import qualified Synthesizer.Plain.Cut as Cut +import qualified Synthesizer.Plain.Signal as Sig + +import qualified Data.EventList.Relative.TimeBody as EventList + +-- import qualified Algebra.RealRing as RealRing +-- import qualified Algebra.Ring as Ring +-- import qualified Algebra.Additive as Additive + +import qualified Number.NonNegative as NonNeg + +import Test.QuickCheck (quickCheck, ) +import Test.Utility (equalList, ) + +import NumericPrelude.Numeric +import NumericPrelude.Base +import Prelude () + + +arrange :: NonNeg.Int -> EventList.T NonNeg.Int (Sig.T Int) -> Bool +arrange nnChunkSize evs = + let chunkSize = SigSt.chunkSize $ 1 + NonNeg.toNumber nnChunkSize + sevs = EventList.mapBody (SigSt.fromList chunkSize) evs + in equalList $ + SigSt.fromList chunkSize (Cut.arrange evs) : + CutSt.arrangeAdaptive chunkSize sevs : + CutSt.arrangeList chunkSize sevs : + CutSt.arrangeEquidist chunkSize sevs : + [] + + +tests :: [(String, IO ())] +tests = + ("arrange", quickCheck arrange) : + [] diff --git a/test/Test/Utility.hs b/test/Test/Utility.hs new file mode 100644 --- /dev/null +++ b/test/Test/Utility.hs @@ -0,0 +1,69 @@ +{-# LANGUAGE NoImplicitPrelude #-} +module Test.Utility where + +import Test.QuickCheck (Arbitrary(arbitrary)) + +import qualified Number.Complex as Complex + +import qualified Algebra.RealRing as RealRing + +import qualified Data.List.HT as ListHT +import qualified Data.Char as Char + +import NumericPrelude.Base +import NumericPrelude.Numeric + + +equalList :: Eq a => [a] -> Bool +equalList xs = + and (ListHT.mapAdjacent (==) xs) + + +approxEqual :: (RealRing.C a) => a -> a -> a -> Bool +approxEqual eps x y = + 2 * abs (x-y) <= eps * (abs x + abs y) + +approxEqualAbs :: (RealRing.C a) => a -> a -> a -> Bool +approxEqualAbs eps x y = + abs (x-y) <= eps + +approxEqualListRel :: (RealRing.C a) => a -> [a] -> Bool +approxEqualListRel eps xs = + let n = fromIntegral $ length xs + in approxEqualListAbs (eps * n * sum (map abs xs)) xs + +approxEqualListAbs :: (RealRing.C a) => a -> [a] -> Bool +approxEqualListAbs eps xs = + let n = fromIntegral $ length xs + s = sum xs + in sum (map (\x -> abs (n*x-s)) xs) <= eps + + +approxEqualComplex :: + (RealRing.C a) => + a -> Complex.T a -> Complex.T a -> Bool +approxEqualComplex eps x y = + 2 * Complex.magnitudeSqr (x-y) + <= eps^2 * (Complex.magnitudeSqr x + Complex.magnitudeSqr y) + +approxEqualComplexAbs :: + (RealRing.C a) => + a -> Complex.T a -> Complex.T a -> Bool +approxEqualComplexAbs eps x y = + Complex.magnitudeSqr (x-y) <= eps^2 + + +-- see event-list + +newtype ArbChar = ArbChar Char + deriving (Eq, Ord) + +instance Show ArbChar where + showsPrec n (ArbChar c) = showsPrec n c + +instance Arbitrary ArbChar where + arbitrary = fmap (ArbChar . Char.chr . (32+) . flip mod 96) arbitrary + +unpackArbString :: [ArbChar] -> String +unpackArbString = + map (\(ArbChar c) -> c)