module Combinatorics.Partitions (
   pentagonalPowerSeries,
   numPartitions,
   partitionsInc,
   partitionsDec,
   allPartitionsInc,

   propInfProdLinearFactors,
   propPentagonalPowerSeries,
   propPentagonalsDifP,
   propPentagonalsDifN,
   ) where

import qualified PowerSeries as PS
import Data.Eq.HT (equating)


{- $setup
>>> import qualified Combinatorics.Partitions as Parts
>>> import qualified Test.QuickCheck as QC
>>> import Data.List (genericLength)
>>> import Data.Eq.HT (equating)
-}


{-
  a(n) denotes the number in how many ways n can be presented as a sum of
  positive integers:
  a(n) n
    1  1 : 1
    2  2 : 2, 1+1
    3  3 : 3, 2+1, 1+1+1
    5  4 : 4, 3+1, 2+2, 2+1+1, 1+1+1+1
    7  5 : 5, 4+1, 3+2, 3+1+1, 2+2+1, 2+1+1+1, 1+1+1+1+1

  Number of partitions: http://oeis.org/A000041
  Pentagonal numbers: http://oeis.org/A001318
-}

{- |
Pentagonal numbers are used to simplify the infinite product
\\prod_{i>0} (1-t^i)
It is known that the coefficients of the power series
are exclusively -1, 0 or 1.
The following is a very simple but inefficient implementation,
because of many multiplications with zero.
-}
prodLinearFactors :: Int -> PS.T Integer
prodLinearFactors :: Int -> [Integer]
prodLinearFactors Int
n =
   forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl forall a. Num a => [a] -> [a] -> [a]
PS.mul [Integer
1] forall a b. (a -> b) -> a -> b
$ forall a. Int -> [a] -> [a]
take Int
n forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map (Integer
1forall a. a -> [a] -> [a]
:) forall a b. (a -> b) -> a -> b
$ forall a. (a -> a) -> a -> [a]
iterate (Integer
0forall a. a -> [a] -> [a]
:) [-Integer
1]

infProdLinearFactors :: PS.T Integer
infProdLinearFactors :: [Integer]
infProdLinearFactors =
   forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith forall a. [a] -> Int -> a
(!!)
      (forall b a. (b -> a -> b) -> b -> [a] -> [b]
scanl (\[Integer]
prod Int
i -> [Integer] -> Int -> [Integer] -> [Integer]
delayedSub [Integer]
prod Int
i [Integer]
prod) [Integer
1] [Int
1..])
      [Int
0..]

{- |
prop> QC.forAll (QC.choose (0,100)) Parts.propInfProdLinearFactors
-}
propInfProdLinearFactors :: Int -> Bool
propInfProdLinearFactors :: Int -> Bool
propInfProdLinearFactors Int
n =
   forall (t :: * -> *). Foldable t => t Bool -> Bool
and forall a b. (a -> b) -> a -> b
$
   forall a. Int -> [a] -> [a]
take (Int
nforall a. Num a => a -> a -> a
+Int
1) forall a b. (a -> b) -> a -> b
$
   forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith forall a. Eq a => a -> a -> Bool
(==)
      [Integer]
infProdLinearFactors
      (Int -> [Integer]
prodLinearFactors Int
n)


pentagonalsP, pentagonalsN,
  pentagonalsDifP, pentagonalsDifN :: [Int]

pentagonalsP :: [Int]
pentagonalsP = forall a b. (a -> b) -> [a] -> [b]
map (\Int
n -> forall a. Integral a => a -> a -> a
div (Int
nforall a. Num a => a -> a -> a
*(Int
3forall a. Num a => a -> a -> a
*Int
nforall a. Num a => a -> a -> a
-Int
1)) Int
2) [Int
0..]
pentagonalsN :: [Int]
pentagonalsN = forall a b. (a -> b) -> [a] -> [b]
map (\Int
n -> forall a. Integral a => a -> a -> a
div (Int
nforall a. Num a => a -> a -> a
*(Int
3forall a. Num a => a -> a -> a
*Int
nforall a. Num a => a -> a -> a
+Int
1)) Int
2) [Int
0..]

{-
  (n+1)*(3*n+2) - n*(3*n-1) = 6*n+2
  (n+1)*(3*n+4) - n*(3*n+1) = 6*n+4
-}
pentagonalsDifP :: [Int]
pentagonalsDifP = forall a b. (a -> b) -> [a] -> [b]
map (\Int
n -> Int
3forall a. Num a => a -> a -> a
*Int
nforall a. Num a => a -> a -> a
+Int
1) [Int
0..]
pentagonalsDifN :: [Int]
pentagonalsDifN = forall a b. (a -> b) -> [a] -> [b]
map (\Int
n -> Int
3forall a. Num a => a -> a -> a
*Int
nforall a. Num a => a -> a -> a
+Int
2) [Int
0..]

{- |
prop> Parts.propPentagonalsDifP 10000
-}
propPentagonalsDifP :: Int -> Bool
propPentagonalsDifP :: Int -> Bool
propPentagonalsDifP Int
n =
   forall b a. Eq b => (a -> b) -> a -> a -> Bool
equating (forall a. Int -> [a] -> [a]
take Int
n)
      [Int]
pentagonalsDifP (forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith (-) (forall a. [a] -> [a]
tail [Int]
pentagonalsP) [Int]
pentagonalsP)

{- |
prop> Parts.propPentagonalsDifN 10000
-}
propPentagonalsDifN :: Int -> Bool
propPentagonalsDifN :: Int -> Bool
propPentagonalsDifN Int
n =
   forall b a. Eq b => (a -> b) -> a -> a -> Bool
equating (forall a. Int -> [a] -> [a]
take Int
n)
      [Int]
pentagonalsDifN (forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith (-) (forall a. [a] -> [a]
tail [Int]
pentagonalsN) [Int]
pentagonalsN)

{-
  delay y by del and subtract it from x
-}
delayedSub :: [Integer] -> Int -> [Integer] -> [Integer]
delayedSub :: [Integer] -> Int -> [Integer] -> [Integer]
delayedSub [Integer]
x Int
del [Integer]
y =
   let ([Integer]
a,[Integer]
b) = forall a. Int -> [a] -> ([a], [a])
splitAt Int
del [Integer]
x
   in  [Integer]
a forall a. [a] -> [a] -> [a]
++ forall a. Num a => [a] -> [a] -> [a]
PS.sub [Integer]
b [Integer]
y

{-
  p00 p01 p02 p03 p04 p05 p06 p07 p08 p09 p10 p11 p12 p13 p14 p15 p16 p17
 -    p00 p01 p02 p03 p04 p05 p06 p07 p08 p09 p10 p11 p12 p13 p14 p15 p16
 +                    p00 p01 p02 p03 p04 p05 p06 p07 p08 p09 p10 p11 p12
 -                                                p00 p01 p02 p03 p04 p05
  ...
 -        p00 p01 p02 p03 p04 p05 p06 p07 p08 p09 p10 p11 p12 p13 p14 p15
 +                            p00 p01 p02 p03 p04 p05 p06 p07 p08 p09 p10
 -                                                            p00 p01 p02
  ...
-}
numPartitions :: [Integer]
numPartitions :: [Integer]
numPartitions =
   let accu :: [Int] -> [Integer]
accu = forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr ([Integer] -> Int -> [Integer] -> [Integer]
delayedSub [Integer]
numPartitions) (forall a. HasCallStack => [Char] -> a
error [Char]
"never evaluated")
       ps :: [Integer]
ps   = [Int] -> [Integer]
accu (forall a. [a] -> [a]
tail [Int]
pentagonalsDifP)
       ns :: [Integer]
ns   = [Int] -> [Integer]
accu (forall a. [a] -> [a]
tail [Int]
pentagonalsDifN)
   in  Integer
1 forall a. a -> [a] -> [a]
: forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith forall a. Num a => a -> a -> a
(+) [Integer]
ps (Integer
0forall a. a -> [a] -> [a]
:[Integer]
ns)

{- |
This is a very efficient implementation of 'prodLinearFactors'.
-}
pentagonalPowerSeries :: [Integer]
pentagonalPowerSeries :: [Integer]
pentagonalPowerSeries =
   let make :: [Int] -> [Integer]
make = forall (t :: * -> *) a. Foldable t => t [a] -> [a]
concat forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith (\Integer
s Int
n -> Integer
s forall a. a -> [a] -> [a]
: forall a. Int -> a -> [a]
replicate (Int
nforall a. Num a => a -> a -> a
-Int
1) Integer
0) (forall a. [a] -> [a]
cycle [Integer
1,-Integer
1])
   in  forall a b c. (a -> b -> c) -> b -> a -> c
flip forall a. Num a => [a] -> [a] -> [a]
PS.sub [Integer
1] forall a b. (a -> b) -> a -> b
$
       forall a. Num a => [a] -> [a] -> [a]
PS.add
          ([Int] -> [Integer]
make [Int]
pentagonalsDifP)
          ([Int] -> [Integer]
make [Int]
pentagonalsDifN)

{- |
prop> Parts.propPentagonalPowerSeries 1000
-}
propPentagonalPowerSeries :: Int -> Bool
propPentagonalPowerSeries :: Int -> Bool
propPentagonalPowerSeries Int
n =
   forall b a. Eq b => (a -> b) -> a -> a -> Bool
equating (forall a. Int -> [a] -> [a]
take Int
n) [Integer]
infProdLinearFactors [Integer]
pentagonalPowerSeries



{- | Give all partitions of the natural number n
     with summands which are at least k.
     Not quite correct for k>n. -}
partitionsInc :: (Integral a) => a -> a -> [[a]]
partitionsInc :: forall a. Integral a => a -> a -> [[a]]
partitionsInc a
k a
n =
   forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\a
y -> forall a b. (a -> b) -> [a] -> [b]
map (a
yforall a. a -> [a] -> [a]
:) (forall a. Integral a => a -> a -> [[a]]
partitionsInc a
y (a
nforall a. Num a => a -> a -> a
-a
y))) [a
k .. forall a. Integral a => a -> a -> a
div a
n a
2] forall a. [a] -> [a] -> [a]
++ [[a
n]]

partitionsDec :: (Integral a) => a -> a -> [[a]]
partitionsDec :: forall a. Integral a => a -> a -> [[a]]
partitionsDec a
0 a
0 = [forall a. a -> [a]
repeat a
0]
partitionsDec a
_ a
0 = []
partitionsDec a
k a
n =
   (if a
kforall a. Ord a => a -> a -> Bool
>=a
n then [[a
n]] else []) forall a. [a] -> [a] -> [a]
++
      forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\a
y -> forall a b. (a -> b) -> [a] -> [b]
map (a
yforall a. a -> [a] -> [a]
:) (forall a. Integral a => a -> a -> [[a]]
partitionsDec a
y (a
nforall a. Num a => a -> a -> a
-a
y)))
                (forall a. (a -> Bool) -> [a] -> [a]
takeWhile (forall a. Ord a => a -> a -> Bool
>a
0) (forall a. (a -> a) -> a -> [a]
iterate forall a. Enum a => a -> a
pred (forall a. Ord a => a -> a -> a
min a
n a
k)))

_partitionsInc :: (Integral a) => a -> a -> [[a]]
_partitionsInc :: forall a. Integral a => a -> a -> [[a]]
_partitionsInc a
k a
n =
   if a
kforall a. Ord a => a -> a -> Bool
>a
n
     then []
     else forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\a
y -> forall a b. (a -> b) -> [a] -> [b]
map (a
yforall a. a -> [a] -> [a]
:) (forall a. Integral a => a -> a -> [[a]]
_partitionsInc a
y (a
nforall a. Num a => a -> a -> a
-a
y))) [a
k..(a
nforall a. Num a => a -> a -> a
-a
1)]
            forall a. [a] -> [a] -> [a]
++ [[a
n]]

{- |
type Int is needed because of list node indexing

prop> QC.forAll (QC.choose (1,10)) $ \k -> QC.forAll (QC.choose (0,50)) $ \n -> Parts.partitionsInc k n == Parts.allPartitionsInc !! k !! n
prop> equating (take 30) (map genericLength (Parts.allPartitionsInc !! 1)) Parts.numPartitions
-}
allPartitionsInc :: [[[[Int]]]]
allPartitionsInc :: [[[[Int]]]]
allPartitionsInc =
   let part :: Int -> Int -> [[Int]]
       part :: Int -> Int -> [[Int]]
part Int
k Int
n = forall (t :: * -> *) a b. Foldable t => (a -> [b]) -> t a -> [b]
concatMap (\Int
y -> forall a b. (a -> b) -> [a] -> [b]
map (Int
yforall a. a -> [a] -> [a]
:) ([[[[Int]]]]
xs forall a. [a] -> Int -> a
!! Int
y forall a. [a] -> Int -> a
!! (Int
nforall a. Num a => a -> a -> a
-Int
y)))
                            [Int
k .. forall a. Integral a => a -> a -> a
div Int
n Int
2]
                      forall a. [a] -> [a] -> [a]
++ [[Int
n]]
       xs :: [[[[Int]]]]
xs = forall a. a -> [a]
repeat [[]] forall a. a -> [a] -> [a]
: forall a b. (a -> b) -> [a] -> [b]
map (\Int
k -> forall a b. (a -> b) -> [a] -> [b]
map (Int -> Int -> [[Int]]
part Int
k) [Int
0..]) [Int
1..]
   in  [[[[Int]]]]
xs