hubFS: THE place for F#

I am trying to do the following:

I have two sequences.  The first is a sequence of all prime numbers. 

The second is a sequence of all numbers (not necessarily prime) that satisfy some property which depends on the primes from the first sequence.  Specifically, it's a sequence of all numbers that do NOT have a prime factor that is congruent to anything other than 1 mod 4.  So, for example, 25 is in the second sequence because 25=5*5, and 5 is congruent to 1 mod 4.  11 is not in the sequence, because 11 is congruent to 3 mod 4, 15 is not in the sequence, because 15=5*3, and even though 5 is congruent to 1 mod 4, 3 is not, and 5525 is in the sequence, because it is equal to 5*5*13*17, all of which are congruent to 1 mod 4. 


I actually don't care about the list of primes, they are only there in order to help me determine if a number is in the second sequence.  Furthermore, I do not know in advance how many values I will end up needing from the second sequence, and on top of that I definitely don't want to test the same number for primality twice, and I will need to make multiple successively deepening passes over the second sequence, so I don't want to ever compute one of those values more than once either.

Is there an elegant solution to this?  I still have trouble thinking "lazily" but it seems like a good candidate.

Also, I'm not too concerned about the algorithms themselves, mostly just how to express the relationship between the two in a way that minimizes recomputation.

Thanks

This does sound like a good candidate for using LazyList.

So long as you have an algorithm for 'compute the next value in the series' (which is what you need to make a LazyList), you should be good to go, then can just treat the list as though it were infinitely populated, but only computes as needed and caches results.

Do not use 'seq', it is pretty useless for this.

so let's say I have

let rec next_prime prev =
    let is_prime n =
        let limit = System.Math.Sqrt(float(n)) |> int
        List.for_all (fun k -> (gcd k n) = 1) [1..limit]
   let this = prev + (2 - prev%2)
   if (is_prime this) then this else next_prime this

Not the best algorithm in the world, but it's short at least for the purposes of illustration.

How can I make a LazyList out of this?  I guess the trouble comes from the fact that when I'm specifying the function for the LazyList via consf, it takes a unit and returns a new LazyList.  How do I access "the last already computed value"?  Also, I assume that any previously computed value will always remain computed correct, so that iterating over it a second time will not force a recomputation?

Here's some code to get you started.

#light

let NextEven prev =
    printfn "NextEven %d" prev
    prev + 2

let rec EvensStartingFrom n = LazyList.consf n (fun () -> EvensStartingFrom (NextEven n))

let evens = EvensStartingFrom 2

let rec WalkPortion n ll =
    if n > 0 then
        match ll with
        | LazyList.Cons(h,t) -> printfn "walked %d" h; WalkPortion (n-1) t
        | LazyList.Nil -> printfn "reached end"

WalkPortion 3 evens
WalkPortion 6 evens

Ok based on that I came up with the following:

let rec gcd a b = if (b=0) then a else gcd b (a%b)
let next_prime prev =
    printfn "Calculating prime after %d" prev
    let rec next_prime prev =
        let is_prime n =
            let limit = System.Math.Sqrt(float(n)) |> int
            List.for_all (fun k -> (gcd k n) = 1) [1..limit]
        let this = if (prev%2<>1) then (prev+1) else (prev+2)
        if (is_prime this) then this else next_prime this
    next_prime prev
let next_k prev primes =
    printfn "Calculating k after %d" prev
    let rec next_k prev primes =
        let this = if (prev%2<>1) then (prev+1) else (prev+2)
        let rec check_all_primes prime_list =
            match prime_list with
            | LazyList.Cons(p,next) ->
                if (p=this) then (p%4=1)
                else if (p>this) then true
                else if (this%p=0 && p%4<>1) then false
                else check_all_primes next
            | LazyList.Nil -> failwith "Prime sequence ends!"
        if (check_all_primes primes) then this
        else next_k this primes
    next_k prev primes
let rec PrimesStartingFrom p =
    LazyList.consf p (fun () -> PrimesStartingFrom (next_prime p))
let KsStartingFrom k =
    let primes = PrimesStartingFrom 2
    let rec KsStartingFrom k =
        LazyList.consf k (fun () -> KsStartingFrom (next_k k primes))
    KsStartingFrom k
let primes = PrimesStartingFrom 2
let ks = KsStartingFrom 1

Is there a way to improve this, or is this pretty much right on?  It does output the right sequence of printfs when I walk the list of k's, and never recomputes anything.

Thanks!

I'm getting the hang of it :P

I decided to go ahead and optimize the algorithm, since when testing n for primality it's really stupid to check every single number <= sqrt(n) for divisibility.  So the idea was to only check primes  less than or equal to sqrt(n).  I struggled over this for a bit because I was dealing with the problem of referencing the list of primes from itself.  I came up with a solution but I get tons of warning about recursive references being checked at runtime for soundness.  What does this mean, and do I need to worry about it?  Also, I feel like the code could be made more elegant.  I'm not sure I like the idea of separating out the is_prime, gen_next, and next_prime functions, it seems like they could somehow be declared inline inside the scope of the definition of the original data structure.  I did it this way though because like I said, the algorithms needed to know about smaller primes in order to determine if larger numbers are prime, and therefore also to specify the function that generates the tail.  So this is what I have:

let is_prime n primes =
    let limit = System.Math.Sqrt(float(n)) |> int
   
    let rec is_prime primes =
        match primes with
        | LazyList.Cons(p,tail) when (p>limit) -> (*printfn "p=%d and limit=%d.  Returning true" p limit;*) true
        | LazyList.Cons(p,tail) ->
            //printfn "Checking if %d is divisible by prime %d" n p
            if (n=p) then (*printfn "It's equal to p, so it's prime";*) true
            else if (n%p=0) then (*printfn "It's divisible, so not prime";*) false
            else is_prime tail
    is_prime primes
   
let next_prime_l prev primes =
    let next = if (prev%2=0) then prev+1 else prev+2
   
    let rec next_prime_l nextcand =
        //printfn "Checking if %d is prime" nextcand
        if (is_prime nextcand primes) then
            //printfn "It is"
            nextcand
        else
            //printfn "It's not, moving on..."
            next_prime_l (nextcand+2)
    next_prime_l next
   
let gen_next primes =
    let rec gen_next tail =
        match tail with
        | LazyList.Cons(hd,tl) ->
            let next = next_prime_l hd primes
            LazyList.consf next (fun () -> gen_next tl)
    gen_next primes
let rec all_primes = LazyList.consf 2 ( fun() -> gen_next all_primes )


Is there a way to do this any better, (for that matter is it even correct)? 

Here's a tighter version, I'm sure you can improve it more.

#light

let rec is_prime n =
    let limit = System.Math.Sqrt(float(n)) |> int
    let rec is_prime primes =
        match primes with
        | LazyList.Cons(p,tail) when (p>limit) ->
            (*printfn "p=%d and limit=%d.  Returning true" p limit;*) true
        | LazyList.Cons(p,tail) ->
            //printfn "Checking if %d is divisible by prime %d" n p
            if (n=p) then (*printfn "It's equal to p, so it's prime";*) true
            else if (n%p=0) then (*printfn "It's divisible, so not prime";*) false
            else is_prime tail
        | _ -> failwith "impossible, primes are infinite list"
    is_prime primes
   
and next_prime_after n =
    assert(n%2=1)
    let candidate = n + 2
    //printfn "Checking if %d is prime" candidate
    if (is_prime candidate) then
        //printfn "It is"
        candidate
    else
        //printfn "It's not, moving on..."
        next_prime_after candidate
   
and all_primes_after n =
    let next = next_prime_after n
    LazyList.consf next (fun () -> all_primes_after next)

and primes = LazyList.consf 2 ( fun() -> LazyList.consf 3 ( fun() -> all_primes_after 3 ) )

primes |> LazyList.take 10 |> Seq.to_list |> printfn "%A"

The warning about initialization soundness is just telling you that it will check for bad programs such as this one:

let Call f = f ()

let rec Kaboom = Call (fun () -> 1 + Kaboom)

which tries to evaluate Kaboom inside the definition of Kaboom.  'primes' is similar, except we don't evaluate primes now, we just store it away in a thunk inside the lazy list tail, so it will only be evaluated after the initial value of 'primes' has been set.