Extensible prime generator: Difference between revisions
→OCaml: add
Antoni Gual (talk | contribs) |
(→OCaml: add) |
||
Line 3,424:
As noted in that section, using an extensible generator isn't the best choice for huge ranges as much higher efficiency can be obtained by using functions that deal directly with the composite culled packed-bit array representations directly as the `countPrimesTo` function does. This makes further improvements to the culling efficiency pointless, as even if one were to reduce the culling speed to zero, it still would take several seconds per range of a billion to enumerate the found primes as a generator.
=={{header|
;2-3-5-7-wheel postponed incremental sieve
<syntaxhighlight lang="ocaml">module IntMap = Map.Make(Int)
let rec steps =
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 :: 2 :: steps
let not_in_wheel =
let scan i =
let rec loop n w = n < 223 && (i = n mod 210
|| match w with [] -> assert false | d :: w' -> loop (n + d) w')
in not (loop 13 steps)
in Array.init 210 scan
let seq_primes =
let rec calc ms m p2 =
if not_in_wheel.(m mod 210) || IntMap.mem m ms
then calc ms (m + p2) p2
else IntMap.add m p2 ms
in
let rec next c p pp ps whl ms () =
match whl with
| [] -> assert false
| d :: w -> match IntMap.min_binding_opt ms with
| Some (m, p2) when c = m ->
next (c + d) p pp ps w (calc (IntMap.remove m ms) (m + p2) p2) ()
| _ when c < pp -> Seq.Cons (c, next (c + d) p pp ps w ms)
| _ -> match ps () with
| Seq.Cons (p', ps') -> let p2' = p + p in
next (c + d) p' (p' * p') ps' w (calc ms (pp + p2') p2') ()
| _ -> assert false
in
let rec ps () = next 13 11 121 ps steps IntMap.empty () in
Seq.cons 2 (Seq.cons 3 (Seq.cons 5 (Seq.cons 7 (Seq.cons 11 ps))))</syntaxhighlight>
;Test code
<syntaxhighlight lang="ocaml">let seq_show sq =
print_newline (Seq.iter (Printf.printf " %u") sq)
let () =
seq_primes |> Seq.take 20 |> seq_show;
seq_primes |> Seq.drop_while ((>) 100) |> Seq.take_while ((>) 150) |> seq_show;
seq_primes |> Seq.drop_while ((>) 7700) |> Seq.take_while ((>) 8000)
|> Seq.length |> Printf.printf " %u primes\n";
seq_primes |> Seq.drop 9999 |> Seq.take 1 |> seq_show</syntaxhighlight>
{{out}}
<pre>
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71
101 103 107 109 113 127 131 137 139 149
30 primes
104729
</pre>
=={{header|PARI/GP}}==
PARI includes a nice prime generator which is quite extensible:
<syntaxhighlight lang="c">void
| |||