fork download 
  1. //Incremental Sieve of Erotasthenes based on the DotNet hash based mutable Dictionary class
  2. //Similar to http://d...content-available-to-author-only...h.net/CommentView,guid,82535152-f7e2-4051-afd6-2330aea0bdf9.aspx#commentstart
  3. //but with opimizations of wheel factorization, delayed adding base primes, and a more efficient use of sequences
  4.  
  5. //#time
  6. [<RequireQualifiedAccess>]
  7. module Dict =
  8.   let empty = null
  9.   let inline tryFind k (m:System.Collections.Generic.Dictionary<_,_>) =
  10.     if m = null then None else let (found,v) = m.TryGetValue k in if found then Some v else None
  11.   let inline add k v (m:System.Collections.Generic.Dictionary<_,_>) =
  12.     let nm = if m = null then new System.Collections.Generic.Dictionary<_,_>() else m in nm.[k] <- v; nm
  13.   let inline remove k (m:System.Collections.Generic.Dictionary<_,_>) = m.Remove k |> ignore; m
  14.  
  15. type cullstate = struct val p:uint32 val pi:int new(p,pi) = { p=p;pi=pi } end
  16. type CIS<'T> = class val v:'T val cont:unit->CIS<'T> new(v,cont) = { v=v;cont=cont } end
  17. let primesDWSE() =
  18.   let WHLPRMS = [| 2u;3u;5u;7u |] in let FSTPRM = 11u in let WHLCRC = int (WHLPRMS |> Seq.fold (*) 1u) >>> 1
  19.   let WHLLMT = int (WHLPRMS |> Seq.fold (fun o n->o*(n-1u)) 1u) - 1
  20.   let WHLPTRN =
  21.     let wp = Array.zeroCreate (WHLLMT+1)
  22.     let gaps (a:int[]) = let rec gap i acc = if a.[i]=0 then gap (i+1) (acc+1uy) else acc
  23.                          {0..WHLCRC-1} |> Seq.fold (fun s i->
  24.                            let ns = if a.[i]<>0 then wp.[s]<-2uy*gap (i+1) 1uy;(s+1) else s in ns) 0 |> ignore
  25.     Array.init (WHLCRC+1) (fun i->if WHLPRMS |> Seq.forall (fun p->(FSTPRM+uint32(i<<<1))%p<>0u)
  26.                                   then 1 else 0) |> gaps;wp
  27.   let inline whladv i = if i < WHLLMT then i + 1 else 0 in let inline advcnd c ci = c + uint32 WHLPTRN.[ci]
  28.   let inline advmltpl c p i = let n = c + uint32 WHLPTRN.[i] * p in if n < c then 0xFFFFFFFFu else n
  29.   let inline reinsert oldcmpst mp (cs:cullstate) =
  30.     let p,pi = cs.p,cs.pi in let cmpst,npi = advmltpl oldcmpst p pi,whladv pi
  31.     match Dict.tryFind cmpst mp with
  32.       | None -> mp |> Dict.add cmpst [cullstate(p,npi)]
  33.       | Some(facts) -> mp |> Dict.add cmpst (cullstate(p,npi)::facts)
  34.   let rec mkprm (n,wi,m,(psd:CIS<_>),q) =
  35.     let nxt,nxti = advcnd n wi,whladv wi
  36.     if nxt < n then (0u,0,(0xFFFFFFFFu,wi,m,psd,q))
  37.     else match Dict.tryFind n m with
  38.            | None -> if n < q then (n,wi,(nxt,nxti,m,psd,q))
  39.                      else let bp,bpi = psd.v in let nc,nci = advmltpl n bp bpi,whladv bpi
  40.                           let nsd = psd.cont() in let np,_ = nsd.v in let sqr = if np>65535u then 0xFFFFFFFFu else np*np
  41.                           mkprm (nxt,nxti,(Dict.add nc [cullstate(bp,nci)] m),nsd,sqr)
  42.            | Some(skips) -> let adjdmap = skips |> List.fold (reinsert n) (m |> Dict.remove n)
  43.                             mkprm (nxt,nxti,adjdmap,psd,q)
  44.   let rec pCIS p pi cont = let ap,api = advcnd p pi,whladv pi
  45.                            CIS((p,pi),fun()->let (np,npi,ncont)=mkprm cont in pCIS np npi ncont)
  46.   let rec baseprimes() = CIS((FSTPRM,0),fun()->let np,npi = advcnd FSTPRM 0,whladv 0
  47.                                                pCIS np npi (advcnd np npi,whladv npi,Dict.empty,baseprimes(),FSTPRM*FSTPRM))
  48.   let inline genseq sd = Seq.unfold (fun (p,_,cont) -> Some(p,mkprm cont)) sd
  49.   seq { yield! WHLPRMS; yield! mkprm (FSTPRM,0,Dict.empty,baseprimes(),(FSTPRM*FSTPRM)) |> genseq }
  50.  
  51. primesDWSE() |> Seq.nth 1000000 |> printfn "%A"
Success #stdin #stdout 3.08s 12864KB
comments (?)
stdin 
Standard input is empty
stdout 
15485867u
https://ideone.com/vjviLS
Optimized DotNet hash table mutable Dictionary class based Sieve of Eratosthenes
language:
F# (mono 4.1)
created:
11 years ago
visibility:
secret

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!