Haskell implementation of some functions inspired by N. J. Cutland's book, ‘Computability. An introduction to recursive function theory’, including the universal function ψ.

URM.hs

module Main () where

import           Prelude hiding ((+), (-), (*), (/), (^), succ, div, min, max)
import qualified Prelude as P
import           Debug.Trace (trace)

optimization = 42


class (Show v) => Val v where
    zero :: v -> Integer
    succ :: v -> Integer
    proj :: Int -> v -> Integer
    len  :: v -> Int

instance Val Integer where
    zero x = 0
    succ x = x P.+ 1
    proj 1 x = x
    proj _ _ = error "wrong proj"
    len _ = 1

instance (Integral n, Show n) => Val [n] where
    zero x       = 0
    succ (x : _) = fromIntegral (x P.+ 1)
    proj p x     = fromIntegral (x !! (p P.- 1))
    len x = length x
    
instance Val () where
    zero ()   = 0
    succ ()   = error "wrong succ"
    proj _ () = error "wrong proj"
    len x = 0
    
compose :: (Val v) => [v -> Integer] -> ([Integer] -> Integer) -> (v -> Integer)
compose lg f inp = 
    let parRis = map (\fx -> fx inp) lg
    in f parRis





primitiveRecursion :: (Val v) => 
    (v -> Integer) -> (v -> Integer -> Integer -> Integer) ->
    (v -> Integer -> Integer)
primitiveRecursion f g = h where
    h x 0   = f x
    h x yp1 = 
        let y = yp1 P.- 1
        in g x y (h x y)

infixl 6 +
(+) = if optimization >= 1
    then (P.+)
    else primitiveRecursion (\x -> x) (\x y z -> succ z)

infixl 7 *
(*) = if optimization >= 1
    then (P.*)
    else primitiveRecursion (\x -> 0) (\x y z -> z + x)

infixr 8 ^
(^) = if optimization >= 1
    then (P.^)
    else primitiveRecursion (\x -> 1) (\x y z -> z * x)

infixl 6 ∸
(∸) = if optimization >= 1
    then fso 
    else fs
    where
        x `fso` y = P.max (x P.- y) 0
        x `fs` 1 = (primitiveRecursion (\() -> 0) (\() y z -> y)) () x
        x `fs` y = (primitiveRecursion (\x  -> x) (\x  y z -> z ∸ 1)) x y

sg  = primitiveRecursion (\() -> 0) (\() y z -> 1) ()
nsg = primitiveRecursion (\() -> 1) (\() y z -> 0) ()

infixl 6 |-|
x |-| y = (x ∸ y) + (y ∸ x)

fact = primitiveRecursion (\() -> 1) (\() y z -> succ(y) * z) ()

min x y = x ∸ (x ∸ y)
max x y = x + (x ∸ y)

rm = if optimization >= 1
    then fro
    else fr 
    where
        0 `fro` y = y
        x `fro` y = P.rem y x
        fr = primitiveRecursion (\x -> 0) (\x y z -> (z + 1) * sg (x |-| (z + 1)))

qt = if optimization >= 1
    then fdo
    else fd 
    where
        0 `fdo` y = 0
        x `fdo` y = P.div y x
        fd = primitiveRecursion (\x -> 0) (\x y z -> z + nsg (x |-| ((rm x y) + 1)))

div x y = nsg (rm x y)


-- Σ(i<l) f(x y)
fΣ :: Val v => Integer -> (v -> Integer -> Integer) -> v -> Integer
fΣ l f x = primitiveRecursion (\x -> 0) (\x y z -> z + f x y) x l

-- Π(i<l) f(x y)
fΠ :: Val v => Integer -> (v -> Integer -> Integer) -> v -> Integer
fΠ l f x = primitiveRecursion (\x -> 1) (\x y z -> z * f x y) x l

-- µ z≤l f(x)
µl :: Val v => Integer -> (v -> Integer -> Integer) -> v -> Integer
µl l f x = 
    fΣ l (\() v ->
        fΠ (v+1) (\() u -> 
            sg (f x u)
        ) ()
    ) ()

-- µ z f(x)
µ :: Val v => (v -> Integer -> Integer) -> v -> Integer
µ f x = head [z | z <- [0..], (f x z) == 0]


d y = fΣ (y+1) (\() x -> nsg (rm x y)) () 

pr = if optimization >= 2
    then fpo
    else fp 
    where
        fp x = nsg ((d x) |-| 2)
        
        x `dv` y = (y `P.mod` x) == 0
        sr = ceiling . sqrt . fromIntegral
        p 0 = False
        p 1 = False
        p 2 = True
        p x = not (2 `dv` x) && null [d | d <- [3, 5..(sr x)], d `dv` x]
        fpo x = if p x then 1 else 0

pᵪ = primitiveRecursion
    (\() -> 0) 
    (\() y z -> µ
        (\() w -> if w > z && (pr w) == 1 then 0 else 1) ()) ()

x ⒳ y = µ (\() z -> div ((pᵪ y) ^ (z+1)) x) ()









-- URM instructions and conversion
π x y = 2^x * (2*y+1) ∸ 1
π1 n  = (n+1) ⒳ 1
π2 n  = qt 2 ((qt (2^(π1 n)) (n+1)) ∸ 1)

ν x y z = 2^x * 3^y * (6*z+1) ∸ 1
ν1 n    = (n+1) ⒳ 1
ν2 n    = (n+1) ⒳ 2
ν3 n    = (((n+1) `P.div` (2 ^ ν1 n)) `P.div` (3 ^ ν2 n)) `P.div` 6

τ :: [Integer] -> Integer
τ l = (p 1 l) ∸ 2 where
    p k [a] = (pᵪ k) ^ (a+1)
    p i (a:l) = ((pᵪ i) ^ a) * p (i+1) l

τl n =
    let x `dv` y = (y `P.mod` x) == 0
        nextp x = head [y | y <- [(x+1)..], pr y == 1]
        pn x = length [y | y <- [2..x], pr y == 1]
        
        fatt x = f x 2 where
            f 1 _ = []
            f x d =
                if d `dv` x
                then d : f (x `P.div` d) d
                else f x (nextp d)
            
    in fromIntegral . pn . P.maximum . fatt $ n+2

τi n i
    | i < τl n  = (n+2) ⒳ i
    | otherwise = ((n+2) ⒳ i) ∸ 1
    
    

data I = Z Integer | S Integer | T Integer Integer | J Integer Integer Integer
    deriving (Show)

β (Z n)     =     4 * (n∸1)
β (S n)     = 1 + 4 * (n∸1)
β (T m n)   = 2 + 4 * (π (m∸1) (n∸1))
β (J m n q) = 3 + 4 * (ν (m∸1) (n∸1) (q∸1))

γ p = τ (map β p)

βi x
    | r == 0 = Z (u+1)
    | r == 1 = S (u+1)
    | r == 2 = T (1 + π1 u) (1 + π2 u)
    | r == 3 = J (1 + ν1 u) (1 + ν2 u) (1 + ν3 u)
    where 
        r = rm 4 x
        u = qt 4 x

γi n = take (τl n) [βi (τi n i) | i <- [1..]]













-- Universal function
zArg q = q+1
sArg q = q+1

tArg1 q = (π1 q) + 1
tArg2 q = (π2 q) + 1

jArg1 q = (ν1 q) + 1
jArg2 q = (ν2 q) + 1
jArg3 q = (ν3 q) + 1

czero  c n = qt ((pᵪ n) ^ (c ⒳ n)) c
csucc  c n = c * (pᵪ n)
ctrasf c m n = (pᵪ n) ^ (c ⒳ m) * (czero c n)

change c i
    | r == 0 = czero  c (zArg q)
    | r == 1 = csucc  c (sArg q)
    | r == 2 = ctrasf c (tArg1 q) (tArg2 q)
    | r == 3 = c
    where
        r = rm 4 i
        q = qt 4 i

nextconf e c t
    | 1 <= t && t <= τl e = change c (τi e t)
    | otherwise = c

is c i t
    | r == 3 && ((c ⒳ (jArg1 q)) == (c ⒳ (jArg2 q)))
                = jArg3 i
    | otherwise = t + 1
    where
        r = rm 4 i
        q = qt 4 i

nextinstr e c t
    | (1 <= t && t <= τl e) && (1 <= isucc && isucc <= (τl e))
                = isucc
    | otherwise = 0
    where isucc = is c (τi e t) t

ck e x 0 = product [(pᵪ (fromIntegral i)) ^ (proj i x) | i <- [1..n]]
    where n = len x
ck e x tp1 
    | True = let t = tp1 P.- 1 in
    nextconf e (ck e x t) (jk e x t)

jk e x 0 = 1
jk e x tp1 = let t = tp1 P.- 1 in
    nextinstr e (ck e x t) (jk e x t)

ψ :: (Val v) => Integer -> v -> Integer
ψ = if optimization >= 3
    then fpo
    else fp 
    where
        fp e x = (ck e x (µ (\x t -> jk e x t) x)) ⒳ 1
        
        fpo e x =
            let lp = τl e
                defconf c t = trace (show (c, t) ++ " " ++ show (map (\i -> c ⒳ i) [1..10])) $ 
                    let nc = nextconf  e c t
                        nt = nextinstr e c t
                    in if 1 <= nt && nt <= lp 
                        then defconf nc nt
                        else nc
            in (defconf (ck e x 0) (jk e x 0)) ⒳ 1
                














addition = [J 3 2 5, S 1, S 3, J 1 1 1]
minusone = [J 1 4 9, S 3, J 1 3 7, S 2, S 3, J 1 1 3, T 2 1]
divtwo   = [J 1 2 6, S 3, S 2, S 2, J 1 1 1, T 3 1]