Infinite subsets of natural numbers

Natural numbers in Haskell

We will show how to use Haskell in order to solve mathematical questions about natural numbers. The problems below are motivated by the Project Euler. The sketches or partial solutions are meant to inspire you to perform further investigation of the problems on your own. You may need to do a substantial optimization of the code, in particular of the underlying data structures, in order to get fast solutions working on inputs proposed in the Project Euler. This is an example how Haskell handles infinite objects:

naturals = [1..]
firstTen = take 10 naturals

main = print ( firstTen )

We decided to represent natural numbers as an infinite list. Haskell allows certain operations on infinite objects, in particular selecting elements from a given list.

Prime numbers

Let us analyze an example. A given natural number n is prime if it is bigger or equal 2 and is not divisible by any natural number except 1 and n. As a preparatory step in our definition of prime numbers we define a function ldf. For two given numbers k and n the function ldf returns the smallest divisor of n which is bigger or equal to k.

ldf :: Integer -> Integer -> Integer
ldf k n | rem n k == 0   = k
        | otherwise      = ldf (k+1) n

main = do
    print ( ldf 2 5 )
    print ( ldf 7 36 )

Now notice, that we can select prime numbers from the naturals using the following procedure:

naturals = [1..]
ldf :: Integer -> Integer -> Integer
ldf k n | rem n k == 0   = k
        | k^2 > n        = n
        | otherwise      = ldf (k+1) n

prime :: Integer -> Bool
prime n | n < 1          = error "this is not a positive number"
        | n == 1         = False
        | otherwise      = ldf 2 n == n
       
primes = filter prime naturals

main = do 
    print ( take 10 primes )

Motivated by Problem 7 from the Project Euler list we may easily figure out the 501st prime number

naturals = [1..]
ldf :: Integer -> Integer -> Integer
ldf k n | rem n k == 0   = k
        | k^2 > n        = n
        | otherwise      = ldf (k+1) n

prime :: Integer -> Bool
prime n | n < 1          = error "this is not a positive number"
        | n == 1         = False
        | otherwise      = ldf 2 n == n

primes = filter prime naturals

main = do 
    print ( primes !! 500 )

Exercise. Modify the code and find 500 number which is not divisible by 2,3,5,7,11.

Collatz sequence

Collatz sequence is an interesting sequence defined in recurential terms. The limit behavior of the sequence is not completely understood. However, we may implement the Collatz sequence in Haskell and ask about its properties such as the longest Collatz sequence starting from a given number n. The is a version of Problem 14 from the Project Euler.

import Data.List

collatz :: Integer -> Maybe (Integer, Integer) 
collatz n | n == 1         = Nothing
          | rem n 2 == 0   = Just ( n, div n 2 )    
          | otherwise      = Just ( n, 3*n+1 )
                 
collatzSeq :: Integer -> [Integer]
collatzSeq n = unfoldr collatz n
 
collatzNumbers = [ length ( collatzSeq n ) | n <- [1..10^3] ]
longestCollatz = maximum collatzNumbers

main = do 
    print ( collatzSeq 15 )
    print ( longestCollatz )
    print ( elemIndex longestCollatz collatzNumbers )

Exercise 2. Implement the Syracuse function and find the longest Syracuse sequence for odd integers below 1000.

Declarative programming

The declarative style of programming implies that some problems find immediate solutions in Haskell. As an example, consider Problem 29 from the Project Euler.

import Data.List

powers = [ a^b | a <- [2..100], b <- [2..100] ]
purePowers = nub powers

main = do 
    print ( length powers )
    print ( length purePowers )

Another example, where the declarative programming allows an immediate solution is Problem 31

import Data.List

change = [ [a,b,c,d,e,f,g,h] | a  <- [0..200],b <- [0..100],c <- [0..50],d <- [0..20],e <- [0..10],f <- [0..4],g <- [0..2],h <- [0..1], a+2*b+5*c+10*d+20*e+50*f+100*g+200*h == 200 ]

main = print ( length change )

A bit more involved problem with prime numbers

Different prime factors, Problem 47

import Data.List

naturals = [1..]

ldf :: Integer -> Integer -> Integer
ldf k n | rem n k == 0   = k
        | k^2 > n        = n
        | otherwise      = ldf (k+1) n

prime :: Integer -> Bool
prime n | n < 1          = error "this is not a positive number"
        | n == 1         = False
        | otherwise      = ldf 2 n == n
       
primes = filter prime naturals

divides :: Integer -> Integer -> Bool
divides n m = rem m n == 0

primeFactors :: Integer -> [Integer]
primeFactors x = unfoldr findFactor x
   where
       findFactor 1 = Nothing
       findFactor x = Just (nextFactor, x `div` nextFactor)
           where nextFactor = head $ filter (\f -> f `divides` x) primes

factors :: Integer -> [Integer]
factors n = nub ( primeFactors n ) 

triples = [ n | n <- [1..10^3], length ( factors n ) == 3, length ( factors (n+1) ) == 3, length ( factors (n+2) ) == 3]

pairs = [ n | n <- [1..10^3], length ( factors n ) == 2, length ( factors (n+1) ) == 2]

main = do
    print ( primeFactors 21 )
    print ( primeFactors 32 )
    print ( factors 32 )
    print ( head triples )
    print ( head pairs )