open System
type SeqDesc<'a> = SeqDesc of 'a * (unit -> SeqDesc<'a>) //' a self referring tuple type
let primesMPWSE =
let whlptrn = [| 2;4;2;4;6;2;6;4;2;4;6;6;2;6;4;2;6;4;6;8;4;2;4;2;
4;8;6;4;6;2;4;6;2;6;6;4;2;4;6;2;6;4;2;4;2;10;2;10 |]
let inline adv i = if i < 47 then i + 1 else 0
let reinsert oldcmpst mp (prime,pi) =
let cmpst = oldcmpst + whlptrn.[pi] * prime
match Map.tryFind cmpst mp with
| None -> mp |> Map.add cmpst [(prime,adv pi)]
| Some(facts) -> mp |> Map.add cmpst ((prime,adv pi)::facts)
let rec mkprimes (n,i) m (SeqDesc((np,npi),nsdf) as psd) q =
let nxt = n + whlptrn.[i]
match Map.tryFind n m with
| None -> if n < q then SeqDesc((n,i),fun() -> mkprimes (nxt,adv i) m psd q)
else let (SeqDesc((nhd,x),ntl) as nsd),nxtcmpst = nsdf(),n + whlptrn.[npi] * np
mkprimes (nxt,adv i) (Map.add nxtcmpst [(np,adv npi)] m) nsd (nhd * nhd)
| Some
(skips
) -> let adjdmap
= skips
|> List.
fold (reinsert n
) (m
|> Map.
remove n
) mkprimes (nxt,adv i) adjdmap psd q
let rec prs = SeqDesc((11,0),fun() -> mkprimes (13,1) Map.empty prs 121 )
let genseq sd = Seq.unfold (fun (SeqDesc((n,i),tailfunc)) -> Some(n,tailfunc())) sd
seq { yield 2; yield 3; yield 5; yield 7; yield! mkprimes (11,0) Map.empty prs 121 |> genseq }
primesMPWSE |> Seq.nth 100000 |> printfn "%A"
b3BlbiBTeXN0ZW0KIAp0eXBlIFNlcURlc2M8J2E+ID0gU2VxRGVzYyBvZiAnYSAqICh1bml0IC0+IFNlcURlc2M8J2E+KSAvLycgYSBzZWxmIHJlZmVycmluZyB0dXBsZSB0eXBlCiAKbGV0IHByaW1lc01QV1NFID0KICBsZXQgd2hscHRybiA9IFt8IDI7NDsyOzQ7NjsyOzY7NDsyOzQ7Njs2OzI7Njs0OzI7Njs0OzY7ODs0OzI7NDsyOwogICAgICAgICAgICAgICAgICAgNDs4OzY7NDs2OzI7NDs2OzI7Njs2OzQ7Mjs0OzY7Mjs2OzQ7Mjs0OzI7MTA7MjsxMCB8XQogIGxldCBpbmxpbmUgYWR2IGkgPSBpZiBpIDwgNDcgdGhlbiBpICsgMSBlbHNlIDAKICBsZXQgcmVpbnNlcnQgb2xkY21wc3QgbXAgKHByaW1lLHBpKSA9CiAgICBsZXQgY21wc3QgPSBvbGRjbXBzdCArIHdobHB0cm4uW3BpXSAqIHByaW1lCiAgICBtYXRjaCBNYXAudHJ5RmluZCBjbXBzdCBtcCB3aXRoCiAgICAgIHwgTm9uZSAtPiBtcCB8PiBNYXAuYWRkIGNtcHN0IFsocHJpbWUsYWR2IHBpKV0KICAgICAgfCBTb21lKGZhY3RzKSAtPiBtcCB8PiBNYXAuYWRkIGNtcHN0ICgocHJpbWUsYWR2IHBpKTo6ZmFjdHMpCiAgbGV0IHJlYyBta3ByaW1lcyAobixpKSBtIChTZXFEZXNjKChucCxucGkpLG5zZGYpIGFzIHBzZCkgcSA9CiAgICBsZXQgbnh0ID0gbiArIHdobHB0cm4uW2ldCiAgICBtYXRjaCBNYXAudHJ5RmluZCBuIG0gd2l0aAogICAgICB8IE5vbmUgLT4gaWYgbiA8IHEgdGhlbiBTZXFEZXNjKChuLGkpLGZ1bigpIC0+IG1rcHJpbWVzIChueHQsYWR2IGkpIG0gcHNkIHEpCiAgICAgICAgICAgICAgICBlbHNlIGxldCAoU2VxRGVzYygobmhkLHgpLG50bCkgYXMgbnNkKSxueHRjbXBzdCA9IG5zZGYoKSxuICsgd2hscHRybi5bbnBpXSAqIG5wCiAgICAgICAgICAgICAgICAgICAgIG1rcHJpbWVzIChueHQsYWR2IGkpIChNYXAuYWRkIG54dGNtcHN0IFsobnAsYWR2IG5waSldIG0pIG5zZCAobmhkICogbmhkKQogICAgICB8IFNvbWUoc2tpcHMpIC0+IGxldCBhZGpkbWFwID0gc2tpcHMgfD4gTGlzdC5mb2xkIChyZWluc2VydCBuKSAobSB8PiBNYXAucmVtb3ZlIG4pCiAgICAgICAgICAgICAgICAgICAgICAgbWtwcmltZXMgKG54dCxhZHYgaSkgYWRqZG1hcCBwc2QgcQogIGxldCByZWMgcHJzID0gU2VxRGVzYygoMTEsMCksZnVuKCkgLT4gbWtwcmltZXMgKDEzLDEpIE1hcC5lbXB0eSBwcnMgMTIxICkKICBsZXQgZ2Vuc2VxIHNkID0gU2VxLnVuZm9sZCAoZnVuIChTZXFEZXNjKChuLGkpLHRhaWxmdW5jKSkgLT4gU29tZShuLHRhaWxmdW5jKCkpKSBzZAogIHNlcSB7IHlpZWxkIDI7IHlpZWxkIDM7IHlpZWxkIDU7IHlpZWxkIDc7IHlpZWxkISBta3ByaW1lcyAoMTEsMCkgTWFwLmVtcHR5IHBycyAxMjEgfD4gZ2Vuc2VxIH0KCnByaW1lc01QV1NFIHw+IFNlcS5udGggMTAwMDAwIHw+IHByaW50Zm4gIiVBIg==