Ordinal.hs

import Data.List

type Natural = Integer -- the natural numbers

-- must be manipulated via the smart constructor `ordinal`
-- to ensure Cantor normal forms.
data Ordinal = Zero | Omega [Ordinal]
  deriving Eq

instance Show Ordinal where
  show Zero = "0"
  show (Omega os) = "(" ++ intercalate " + " [ "ω^" ++ show o  | o <- os ] ++ ")"

height :: Ordinal -> Natural
height Zero = 0
height (Omega os) = 1 + maximum [ height o | o <- os ]

gtLex :: [Ordinal] -> [Ordinal] -> Bool
gtLex [] _ = False
gtLex (_ : _) [] = True
gtLex (o1 : os1) (o2 : os2) = o1 > o2 || (o1 == o2 && gtLex os1 os2)

instance Ord Ordinal where
  Zero > _ = False
  Omega _ > Zero = True
  (Omega os1) > (Omega os2) = gtLex os1 os2
  o1 <= o2 = o1 == o2 || o2 > o1

-- smart constructors
zero :: Ordinal
zero = Zero

one :: Ordinal
one = Omega [Zero]

omega :: [Ordinal] -> Ordinal
omega [] = undefined
omega os = Omega (reverse (sort os))

nat :: Natural -> Ordinal
nat n = nsum (replicate (fromIntegral n) one)

-- natural addition #
nadd :: Ordinal -> Ordinal -> Ordinal
nadd Zero o = o
nadd o Zero = o
nadd (Omega os1) (Omega os2) = omega (os1 ++ os2)

-- natural sum
nsum :: [Ordinal] -> Ordinal
nsum os = foldl nadd Zero os

-- >>> nat 10
-- (ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0)
-- >>> nat 0
-- 0
--

-- >>> omega [ nat 3, nat 10 ]
-- (ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0))
--
-- >>> omega [ zero, omega [one] ]
-- (ω^(ω^(ω^0)) + ω^0)
--
-- >>> (omega [omega [omega [zero]]]) < (omega [omega [zero]])
-- False
--

base' :: Natural -> Natural -> Natural -> [Natural]
base' n b _
  | b < 2 = undefined
  | n < 0 = undefined
base' 0 _ _ = []
base' n b acc = d : (base' (n - d * acc) b (b * acc))
  where d = (n `mod` (b * acc)) `div` acc

-- >>> base' 0 10 1
-- []
-- >>> base' 1 10 1
-- [1]
-- >>> base' 1101 10 1
-- [1,0,1,1]
--

base :: Natural -> Natural -> [Natural]
base n b = reverse (base' n b 1)

-- >>> base 1101 10
-- [1,1,0,1]
--
-- >>> base 1023 2
-- [1,1,1,1,1,1,1,1,1,1]
--

ordinalBase :: Natural -> Natural -> Ordinal
ordinalBase n b
  | b < 2 = undefined
  | n == 0 = Zero
  | otherwise = nsum [ nsum (replicate (fromIntegral d) o) | (k, d) <- zip [0..] (reverse (base n b)), let o = omega [ nat k ] ]

-- >>> ordinalBase 10 2
-- (ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0))
--
-- >>> ordinalBase 25739 10
-- (ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0 + ω^0)
--

substitute :: Ordinal -> Natural -> Natural
substitute Zero _ = 0
substitute (Omega os) b = sum [ b ^ (substitute o b) | o <- os ]

-- >>> substitute (ordinalBase 10 2) 2
-- 10
--
-- >>> substitute (ordinalBase 25739 10) 10
-- 25739
--

weakGoodsteinStep :: Natural -> Natural -> Natural
weakGoodsteinStep n b = substitute (ordinalBase n b) (b + 1) -1

-- >>> weakGoodsteinStep 10 2
-- 29
--
-- >>> weakGoodsteinStep 29 3
-- 65
--
-- >>> weakGoodsteinStep 65 4
-- 125
--

weakGoodstein' :: Natural -> Natural -> [Natural]
weakGoodstein' 0 _ = []
weakGoodstein' n b = next : weakGoodstein' next (b + 1)
  where next = weakGoodsteinStep n b

weakGoodstein :: Natural -> [Natural]
weakGoodstein n = n : weakGoodstein' n 2

-- >>> weakGoodstein 5
-- [5,9,15,17,19,21,23,24,25,26,27,28,29,30,31,31,31,31,31,31,31,31,31,31,31,31,31,31,31,31,31,30,29,28,27,26,25,24,23,22,21,20,19,18,17,16,15,14,13,12,11,10,9,8,7,6,5,4,3,2,1,0]
--
-- >>> weakGoodstein 4
-- [4,8,9,10,11,11,11,11,11,11,11,10,9,8,7,6,5,4,3,2,1,0]
--
-- >>> weakGoodstein 3
-- [3,3,3,2,1,0]
--
-- >>> weakGoodstein 6
-- [6,11,17,25,35,39,43,47,51,55,59,62,65,68,71,74,77,80,83,86,89,92,95,97,99,101,103,105,107,109,111,113,115,117,119,121,123,125,127,129,131,133,135,137,139,141,143,144,145,146,147,148,149,150,151,152,153,154,155,156,157,158,159,160,161,162,163,164,165,166,167,168,169,170,171,172,173,174,175,176,177,178,179,180,181,182,183,184,185,186,187,188,189,190,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,191,190,189,188,187,186,185,184,183,182,181,180,179,178,177,176,175,174,173,172,171,170,169,168,167,166,165,164,163,162,161,160,159,158,157,156,155,154,153,152,151,150,149,148,147,146,145,144,143,142,141,140,139,138,137,136,135,134,133,132,131,130,129,128,127,126,125,124,123,122,121,120,119,118,117,116,115,114,113,112,111,110,109,108,107,106,105,104,103,102,101,100,99,98,97,96,95,94,93,92,91,90,89,88,87,86,85,84,83,82,81,80,79,78,77,76,75,74,73,72,71,70,69,68,67,66,65,64,63,62,61,60,59,58,57,56,55,54,53,52,51,50,49,48,47,46,45,44,43,42,41,40,39,38,37,36,35,34,33,32,31,30,29,28,27,26,25,24,23,22,21,20,19,18,17,16,15,14,13,12,11,10,9,8,7,6,5,4,3,2,1,0]
--

showWeakGoodstein :: Natural -> Int -> String
showWeakGoodstein n k = intercalate "\n>\n" [ show (ordinalBase m b) | (b, m) <- zip [2..] (take k (weakGoodstein n)) ]

-- >>> putStr (showWeakGoodstein 10 5)
-- (ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0))
-- >
-- (ω^(ω^0 + ω^0 + ω^0) + ω^0 + ω^0)
-- >
-- (ω^(ω^0 + ω^0 + ω^0) + ω^0)
-- >
-- (ω^(ω^0 + ω^0 + ω^0))
-- >
-- (ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0 + ω^0 + ω^0)
--


hereditaryOrdinalBase :: Natural -> Natural -> Ordinal
hereditaryOrdinalBase n b
  | b < 2 = undefined
  | n < b = nat n
  | otherwise = nsum [ nsum (replicate (fromIntegral d) o) | (k, d) <- zip [0..] (reverse (base n b)), d /= 0, let o = omega [hereditaryOrdinalBase k b] ] -- o^k . d

-- >>> ordinalBase 10 2
-- (ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0))
--
-- >>> hereditaryOrdinalBase 10 2
-- (ω^(ω^(ω^0) + ω^0) + ω^(ω^0))
--
-- >>> substitute (hereditaryOrdinalBase 10 2) 2
-- 10
--

-- >>> hereditaryOrdinalBase 83 3
-- (ω^(ω^(ω^0) + ω^0) + ω^0 + ω^0)
--



goodsteinStep :: Natural -> Natural -> Natural
goodsteinStep n b = substitute (hereditaryOrdinalBase n b) (b + 1) -1

goodstein' :: Natural -> Natural -> [Natural]
goodstein' 0 _ = []
goodstein' n b = next : goodstein' next (b + 1)
  where next = goodsteinStep n b

goodstein :: Natural -> [Natural]
goodstein n = n : goodstein' n 2

--- >>> take 7 (goodstein 10)
--- [10,83,1025,15625,279935,4215754,84073323]
---

--- >>> hereditaryOrdinalBase 10 2
--- (ω^(ω^(ω^0) + ω^0) + ω^(ω^0))
---
--- >>> hereditaryOrdinalBase 4215754 7
--- (ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0 + ω^0)
--- >>>  hereditaryOrdinalBase 10 2 > hereditaryOrdinalBase 4215754 7
--- True

-- >>> take 3 (goodstein 100)
-- [100,228767924549636,3486030061785075245889246499533519993144635113354022278208125975367658647819122213968487317523394891199408283128486373221976009531106108186072748741581865001]

showGoodstein :: Natural -> Int -> String
showGoodstein n k = intercalate "\n>\n" [ show (hereditaryOrdinalBase m b) | (b, m) <- zip [2..] (take k (goodstein n)) ]

-- >>> putStr (showGoodstein 10 7)
-- (ω^(ω^(ω^0) + ω^0) + ω^(ω^0))
-- >
-- (ω^(ω^(ω^0) + ω^0) + ω^0 + ω^0)
-- >
-- (ω^(ω^(ω^0) + ω^0) + ω^0)
-- >
-- (ω^(ω^(ω^0) + ω^0))
-- >
-- (ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0 + ω^0 + ω^0)
-- >
-- (ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0 + ω^0)
-- >
-- (ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^(ω^0)) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0 + ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^(ω^0) + ω^0 + ω^0 + ω^0)
--