Sieve of Eratosthenes: Difference between revisions

{{trans|Python}}

 

V is_prime = [0B]*2 [+] [1B]*(limit - 1)

L(n) 0 .< Int(limit ^ 0.5 + 1.5)

R enumerate(is_prime).filter((i, prime) -> prime).map((i, prime) -> i)

 

 

{{out}}

=={{header|360 Assembly}}==

For maximum compatibility, this program uses only the basic instruction set.

ERATOST CSECT

USING ERATOST,R12

CRIBLE DC 100000X'01'

YREGS

{{out}}

<pre style="height:20ex">

=={{header|6502 Assembly}}==

If this subroutine is called with the value of <i>n</i> in the accumulator, it will store an array of the primes less than <i>n</i> beginning at address 1000 hex and return the number of primes it has found in the accumulator.

LDA #$00

LDX #$00

JMP COPY

COPIED: TYA ; how many found

 

=={{header|68000 Assembly}}==

Some of the macro code is derived from the examples included with EASy68K.

See 68000 "100 Doors" listing for additional information.

* Title : BitSieve

* Written by : G. A. Tippery

Summary2: DC.B ' prime numbers found.',CR,LF,$00

 

 

 

=={{header|8086 Assembly}}==

 

cpu 8086

bits 16

numbuf: db 13,10,'$'

section .bss

 

{{out}}

4993

4999

</pre>

 

=={{header|AArch64 Assembly}}==

{{works with|as|Raspberry Pi 3B version Buster 64 bits}}

/* ARM assembly AARCH64 Raspberry PI 3B */

/* program cribleEras64.s */

/* for this file see task include a file in language AArch64 assembly */

.include "../includeARM64.inc"

<pre>

Prime : 2

 

=={{header|ABAP}}==

PARAMETERS: p_limit TYPE i OBLIGATORY DEFAULT 100.

 

ENDIF.

ENDLOOP.

 

=={{header|ACL2}}==

(declare (xargs :measure (nfix (- n i))))

(if (zp (- n i))

 

(defun sieve (limit)

 

=={{header|ActionScript}}==

Works with ActionScript 3.0 (this is utilizing the actions panel, not a separated class file)

{

var primes:Array = new Array();

}

var e:Array = eratosthenes(1000);

Output:

{{out}}

=={{header|Ada}}==

 

 

procedure Eratos is

Last: Positive := Positive'Value(Ada.Command_Line.Argument(1));

Prime: array(1 .. Last) of Boolean := (1 => False, others => True);

Cnt: Positive;

begin

if Prime(Base) then

Cnt := Base + Base;

Prime(Cnt) := False;

Cnt := Cnt + Base;

end if;

end loop;

 

{{out}}

=={{header|Agda}}==

 

-- imports

open import Data.Nat as ℕ using (ℕ; suc; zero; _+_; _∸_)

remove _ ys zero = nothing ∷ ys i

remove y ys (suc j) = y ∷ ys j

 

=={{header|Agena}}==

Tested with Agena 2.9.5 Win32

 

# generate and return a sequence containing the primes up to sieveSize

 

# test the sieve proc

{{out}}

<pre>

 

'''Works with:''' ALGOL 60 for OS/360

'INTEGER' 'ARRAY' CANDIDATES(/0..1000/);

'INTEGER' I,J,K;

'END'

'END'

 

=={{header|ALGOL 68}}==

PROC eratosthenes = (INT n)[]BOOL:

(

);

{{out}}

<pre>

=={{header|ALGOL W}}==

=== Standard, non-optimised sieve ===

 

% implements the sieve of Eratosthenes %

end

 

{{out}}

<pre>

 

Alternative version that only stores odd numbers greater than 1 in the sieve.

% implements the sieve of Eratosthenes %

% only odd numbers appear in the sieve, which starts at 3 %

for i := 1 until arrayMax do if s( i ) then writeon( ( i * 2 ) + 1 );

end

{{out}}

Same as the standard version.

 

=={{header|ALGOL-M}}==

BEGIN

 

 

END

{{out}}

<pre>

=== Non-Optimized Version ===

 

b←⍵⍴1

b[⍳2⌊⍵]←0

}

 

 

The required list of prime divisors obtains by recursion (<tt>{⍵/⍳⍴⍵}∇⌈⍵*0.5</tt>).

 

 

b←⍵⍴{∧⌿↑(×/⍵)⍴¨~⍵↑¨1}2 3 5

b[⍳6⌊⍵]←(6⌊⍵)⍴0 0 1 1 0 1

}

 

 

The optimizations are as follows:

=== Examples ===

 

2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

 

 

+/∘sieve¨ 10*⍳10

 

The last expression computes the number of primes < 1e0 1e1 ... 1e9.

=={{header|AppleScript}}==

 

script o

property numberList : {missing value}

end sieveOfEratosthenes

 

 

{{out}}

 

=={{header|ARM Assembly}}==

{{works with|as|Raspberry Pi}}

 

/* ARM assembly Raspberry PI */

/***************************************************/

.include "../affichage.inc"

<pre>

Prime : 2

Prime : 101

</pre>

=={{header|AutoHotkey}}==

{{AutoHotkey case}}Source: [http://www.autohotkey.com/forum/topic44657.html AutoHotkey forum] by Laszlo

 

Sieve(n) { ; Sieve of Eratosthenes => string of 0|1 chars, 1 at position k: k is prime

}

Return S

 

===Alternative Version===

arr := []

loop % n-1

}

return Arr

Arr := Sieve_of_Eratosthenes(n)

loop, % n-1

output .= (Arr[A_Index] ? A_Index : ".") . (!Mod(A_Index, 10) ? "`n" : "`t")

MsgBox % output

{{out}}

<pre>. 2 3 . 5 . 7 . . .

 

=={{header|AutoIt}}==

$M = InputBox("Integer", "Enter biggest Integer")

Global $a[$M], $r[$M], $c = 1

$r[0] = $c - 1

ReDim $r[$c]

 

=={{header|AWK}}==

input from commandline as well as stdin,

and input is checked for valid numbers:

# usage: gawk -v n=101 -f sieve.awk

 

 

END { print "Bye!" }

 

Here is an alternate version that uses an associative array to record composites with a prime dividing it. It can be considered a slow version, as it does not cross out composites until needed. This version assumes enough memory to hold all primes up to ULIMIT. It prints out noncomposites greater than 1.

BEGIN { ULIMIT=100

 

else print ( S[(n+n)] = n )

}

 

==Bash==

{{works with|FreeBASIC}}

{{works with|RapidQ}}

 

INPUT "Enter number to search to: "; limit

FOR n = 2 TO limit

IF flags(n) = 0 THEN PRINT n; ", ";

 

==={{header|Applesoft BASIC}}===

20 DIM FLAGS(LIMIT)

30 FOR N = 2 TO SQR (LIMIT)

100 FOR N = 2 TO LIMIT

110 IF FLAGS(N) = 0 THEN PRINT N;", ";

 

==={{header|Atari BASIC}}===

{{trans|Commodore BASIC}}

Auto-initialization of arrays is not reliable, so we have to do our own. Also, PRINTing with commas doesn't quite format as nicely as one might hope, so we do a little extra work to keep the columns lined up.

110 PRINT "LIMIT";:INPUT LI

120 DIM N(LI):FOR I=0 TO LI:N(I)=1:NEXT I

250 IF C=3 THEN PRINT:C=0

260 NEXT I

{{Out}}

<pre> Ready

==={{header|Commodore BASIC}}===

Since C= BASIC initializes arrays to all zeroes automatically, we avoid needing our own initialization loop by simply letting 0 mean prime and using 1 for composite.

110 INPUT "LIMIT";LI

120 DIM N(LI)

230 NEXT I

240 PRINT

{{Out}}

<pre>

 

==={{header|IS-BASIC}}===

110 LET LIMIT=100

120 NUMERIC T(1 TO LIMIT)

230 FOR I=2 TO LIMIT ! Display the primes

240 IF T(I)=0 THEN PRINT I;

 

==={{header|Locomotive Basic}}===

20 INPUT "Limit";limit

30 DIM f(limit)

100 FOR n=2 TO limit

110 IF f(n)=0 THEN PRINT n;",";

 

==={{header|MSX Basic}}===

6 REM Translated from Rosetta's ZX Spectrum implementation

10 INPUT "Enter number to search to: ";l

100 FOR n=2 TO l

110 IF p(n)=0 THEN PRINT n;", ";

 

==={{header|Sinclair ZX81 BASIC}}===

 

A note on <code>FAST</code> and <code>SLOW</code>: under normal circumstances the CPU spends about 3/4 of its time driving the display and only 1/4 doing everything else. Entering <code>FAST</code> mode blanks the screen (which we do not want to update anyway), resulting in substantially improved performance; we then return to <code>SLOW</code> mode when we have something to print out.

20 FAST

30 DIM N(L)

110 FOR I=2 TO L

120 IF NOT N(I) THEN PRINT I;" ";

 

==={{header|ZX Spectrum Basic}}===

20 DIM p(l)

30 FOR n=2 TO SQR l

100 FOR n=2 TO l

110 IF p(n)=0 THEN PRINT n;", ";

 

==={{header|QL SuperBASIC}}===

{{trans|ZX Spectrum Basic}}

Sets h$ to 1 for higher multiples of 2 via <code>FILL$</code>, later on sets <code>STEP</code> to 2n; replaces Floating Pt array p(z) with string variable h$(z) to sieve out all primes < z=441 (<code>l</code>=21) in under 1K, so that h$ is fillable to its maximum (32766), even on a 48K ZX Spectrum if translated back.

10 INPUT "Enter Stopping Pt for squared factors: ";z

15 LET l=SQR(z)

90 REM Display the primes

100 FOR n=2 TO z: IF h$(n)=0: PRINT n;", ";

 

====2i wheel emulation of Sinclair ZX81 BASIC====

10 INPUT Z

15 LET L=SQR(Z)

120 IF H$(I)="0" THEN PRINT I!

130 NEXT I

 

====2i wheel emulation of Sinclair ZX80 BASIC====

. . . with 2:1 compression (of 16-bit integer variables on ZX80s) such that it obviates having to account for any multiple of 2; one has to input odd upper limits on factors to be squared, <code>L</code> (=21 at most on 1K ZX80s for all primes till 439).

10 INPUT L

15 LET Z=(L+1)*(L- 1)/2

120 IF NOT H(I) THEN PRINT 0^I+1+I*2!

130 NEXT I

 

====Sieve of Sundaram====

just a slight transformation away from the Sieve of Sundaram, as transformed as follows: <code>O</code> is the highest value for an Index of succesive diagonal elements in Sundaram's matrix, for which H(J) also includes the off-diagonal elements in-between, such that

duplicate entries are omitted. Thus, a slightly transformed Sieve of Sundaram is what Eratosthenes' Sieve becomes upon applying all

10 INPUT O

15 LET Z=2*O*O+O*2

40 FOR I=1 TO O

45 LET A=2*I*I+I*2

60 FOR J=A TO Z STEP 1+I*2

70 LET H(J)=1

80 NEXT J

120 IF NOT H(I) THEN PRINT 0^I+1+I*2!

130 NEXT I

 

====Eulerian optimisation====

While slower than the optimised Sieve of Eratosthenes before it, the Sieve of Sundaram above has a compatible compression scheme that's more convenient than the conventional one used beforehand. It is therefore applied below along with Euler's alternative optimisation in a reversed implementation that lacks backward-compatibility to ZX80 BASIC. This program is designed around features & limitations of the QL, yet can be rewritten more efficiently for 1K ZX80s, as they allow integer variables to be parameters of <code>FOR</code> statements (& as their 1K of static RAM is equivalent to L1 cache, even in <code>FAST</code> mode). That's left as an exercise for ZX80 enthusiasts, who for o%=14 should be able to generate all primes < 841, i.e. 3 orders of (base 2) magnitude above the limit for the program listed under Sinclair ZX81 BASIC. In QL SuperBASIC, o% may at most be 127--generating all primes < 65,025 (almost 2x the upper limit for indices & integer variables used to calculate them ~2x faster than for floating point as used in line 30, after which the integer code mimics an assembly algorithm for the QL's 68008.)

10 INPUT "Enter highest value of diagonal index q%: ";o%

15 LET z%=o%*(2+o%*2) : h$=FILL$(" ",z%+o%) : q%=1 : q=q% : m=z% DIV (2*q%+1)

110 GOTO 100

127 DEF FN N%(i)=0^i+1+i*2

 

=={{header|Batch File}}==

@echo off

setlocal ENABLEDELAYEDEXPANSION

if !crible.%%i! EQU 1 echo %%i

)

{{Out}}

<pre style="height:20ex">limit: 100

 

=={{header|BBC BASIC}}==

DIM sieve% limit%

FOR I% = 1 TO limit%

IF sieve%?I% = 0 PRINT I%;

 

=={{header|BCPL}}==

 

manifest $( LIMIT = 1000 $)

$)

wrch('*N')

{{out}}

<pre> 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71

=== Odds-only bit packed array version (64 bit) ===

This sieve also uses an iterator structure to enumerate the primes in the sieve. It's inspired by the golang bit packed sieve that returns a closure as an iterator. However, BCPL does not support closures, so the code uses an iterator object.

GET "libhdr"

 

RESULTIS 0

}

{{Out}}

<pre>

@ ^ p3\" ":<

2 234567890123456789012345678901234567890123456789012345678901234567890123456789

 

=={{header|Bracmat}}==

This solution does not use an array. Instead, numbers themselves are used as variables. The numbers that are not prime are set (to the silly value "nonprime"). Finally all numbers up to the limit are tested for being initialised. The uninitialised (unset) ones must be the primes.

= n j i

. !arg:?n

)

& eratosthenes$100

{{out}}

2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

 

=={{header|C}}==

#include <math.h>

 

}

return sieve;

 

'''Another example:'''

Then, in a loop, fill zeroes into those places where i * j is less than or equal to n (number of primes requested), which means they have multiples!

To understand this better, look at the output of the following example.

#include <malloc.h>

void sieve(int *, int);

}

printf("\n\n");

j:2

Before a[2*2]: 1

i:9

i:10

 

=={{header|C sharp|C#}}==

{{works with|C sharp|C#|2.0+}}

using System.Collections;

using System.Collections.Generic;

}

}

 

 

{{trans|F#}}

 

To show that C# code can be written in somewhat functional paradigms, the following in an implementation of the Richard Bird sieve from the Epilogue of [Melissa E. O'Neill's definitive article](http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf) in Haskell:

using System.Collections;

using System.Collections.Generic;

}

IEnumerator IEnumerable.GetEnumerator() { return (IEnumerator)GetEnumerator(); }

 

 

{{trans|F#}}

 

The above code can easily be converted to "'''odds-only'''" and a infinite tree-like folding scheme with the following minor changes:

using System.Collections;

using System.Collections.Generic;

}

IEnumerator IEnumerable.GetEnumerator() { return (IEnumerator)GetEnumerator(); }

 

The above code runs over ten times faster than the original Richard Bird algorithm.

 

 

{{trans|F#}}

 

First, an implementation of a Min Heap Priority Queue is provided as extracted from the entry at [http://rosettacode.org/wiki/Priority_queue#C.23 RosettaCode], with only the necessary methods duplicated here:

using KeyT = System.UInt32;

using System;

public static MinHeapPQ<V> replaceMin(KeyT k, V v, MinHeapPQ<V> pq) {

pq.rplcmin(k, v); return pq; }

 

The following code implements an improved version of the '''odds-only''' O'Neil algorithm, which provides the improvements of only adding base prime composite number streams to the queue when the sieved number reaches the square of the base prime (saving a huge amount of memory and considerable execution time, including not needing as wide a range of a type for the internal prime numbers) as well as minimizing stream processing using fusion:

using System.Collections;

using System.Collections.Generic;

public IEnumerator<PrimeT> GetEnumerator() { return nmrtr(); }

IEnumerator IEnumerable.GetEnumerator() { return (IEnumerator)GetEnumerator(); }

 

The above code is at least about 2.5 times faster than the Tree Folding version.

 

 

The above code adds quite a bit of overhead in having to provide a version of a Priority Queue for little advantage over a Dictionary (hash table based) version as per the code below:

using System.Collections;

using System.Collections.Generic;

public IEnumerator<PrimeT> GetEnumerator() { return nmrtr(); }

IEnumerator IEnumerable.GetEnumerator() { return (IEnumerator)GetEnumerator(); }

 

The above code runs in about three quarters of the time as the above Priority Queue based version for a range of a million primes which will fall even further behind for increasing ranges due to the Dictionary providing O(1) access times as compared to the O(log n) access times for the Priority Queue; the only slight advantage of the PQ based version is at very small ranges where the constant factor overhead of computing the table hashes becomes greater than the "log n" factor for small "n".

 

 

All of the above unbounded versions are really just an intellectual exercise as with very little extra lines of code above the fastest Dictionary based version, one can have an bit-packed page-segmented array based version as follows:

using System.Collections;

using System.Collections.Generic;

public IEnumerator<PrimeT> GetEnumerator() { return nmrtr(); }

IEnumerator IEnumerable.GetEnumerator() { return (IEnumerator)GetEnumerator(); }

 

The above code is about 25 times faster than the Dictionary version at computing the first about 50 million primes (up to a range of one billion), with the actual enumeration of the result sequence now taking longer than the time it takes to cull the composite number representation bits from the arrays, meaning that it is over 50 times faster at actually sieving the primes. The code owes its speed as compared to a naive "one huge memory array" algorithm to using an array size that is the size of the CPU L1 or L2 caches and using bit-packing to fit more number representations into this limited capacity; in this way RAM memory access times are reduced by a factor of from about four to about 10 (depending on CPU and RAM speed) as compared to those naive implementations, and the minor computational cost of the bit manipulations is compensated by a large factor in total execution time.

 

 

In order to make further gains in speed, custom methods must be used to avoid using iterator sequences. If this is done, then further gains can be made by extreme wheel factorization (up to about another about four times gain in speed) and multi-processing (with another gain in speed proportional to the actual independent CPU cores used).

 

All of the above unbounded code can be tested by the following "main" method (replace the name "PrimesXXX" with the name of the class to be tested):

 

To produce the following output for all tested versions (although some are considerably faster than others):

This implementation follows the standard library pattern of [http://en.cppreference.com/w/cpp/algorithm/iota std::iota]. The start and end iterators are provided for the container. The destination container is used for marking primes and then filled with the primes which are less than the container size. This method requires no memory allocation inside the function.

 

#include <algorithm>

#include <vector>

 

{

std::fill(start, end, 0);

{

{

}

}

{

if (*it == 0)

{

}

}

}

 

{

 

=== Boost ===

 

template<class UnaryFunction>

void primesupto(int limit, UnaryFunction yield)

if (is_prime[n])

yield(n);

 

Full program:

 

$ g++ -I/path/to/boost sieve.cpp -o sieve && sieve 10000000

*/

cout << "limit sum pi(n)\n"

<< limit << " " << sum << " " << count << endl;

 

=={{header|Chapel}}==

{{incorrect|Chapel|Doesn't compile since at least Chapel version 1.20 to 1.24.1.}}

This solution uses nested iterators to create new wheels at run time:

iter sieve(prime):int {

 

yield candidate;

}

for p in sieve(2) {

if p > N {

}

write(" ", p);

 

===Alternate Conventional Bit-Packed Implementation===

compile with the `--fast` option

 

use BitOps;

 

for p in 2 .. n do if cmpsts[p >> 3] & (1: uint(8) << (p & 7)) == 0 then yield p;

 

 

{{out}}

===Alternate Odds-Only Bit-Packed Implementation===

 

use BitOps;

 

if cmpsts[i >> 3] & (1: uint(8) << (i & 7)) == 0 then yield i + i + 3;

 

 

{{out}}

 

{{trans|Python}} [https://rosettacode.org/wiki/Sieve_of_Eratosthenes#Infinite_generator_with_a_faster_algorithm code link]

 

 

config const limit = 100000000;

write("Found ", count, " primes up to ", limit);

writeln(" in ", timer.elapsed(TimeUnits.milliseconds), " milliseconds.");

 

{{out}}

 

<pre>The first 25 primes are: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

 

Time as run using Chapel version 24.1 on an Intel Skylake i5-6500 at 3.6 GHz (turbo, single threaded).

 

 

As an alternate to the use of a built-in library, the following code implements a specialized BasePrimesTable that works similarly to the way the Python associative arrays work as to hashing algorithm used (no hashing, as the hash values for integers are just themselves) and something similar to the Python method of handling hash table collisions is used:

 

Compile with the `--fast` compiler command line option

 

config const limit = 100000000;

write("Found ", count, " primes up to ", limit);

writeln(" in ", timer.elapsed(TimeUnits.milliseconds), " milliseconds.");

 

{{out}}

 

 

This last code is quite usable up to a hundred million (as here) or even a billion in a little over ten times the time, but is still slower than the very simple odds-only monolithic array version and is also more complex, although it uses less memory (only for the hash table for the base primes of about eight Kilobytes for sieving to a billion compared to over 60 Megabytes for the monolithic odds-only simple version).

 

===Functional Tree Folding Odds-Only Version===

{{trans|Nim}} [https://rosettacode.org/wiki/Sieve_of_Eratosthenes#Nim_Unbounded_Versions code link]

{{works with|Chapel|1.22|- compile with the --fast compiler command line flag for full optimization}}

 

type Prime = uint(32);

for p in primes() { if p > limit then break; cnt += 1; }

timer.stop(); write("\nFound ", cnt, " primes up to ", limit);

{{out}}

<pre>The first 25 primes are: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

To take advantage of the features that make Chapel shine, we need to use it to do some parallel computations, and to efficiently do that for the Sieve of Eratosthenes, we need to divide the work into page segments where we can assign each largish segment to a separate thread/task; this also improves the speed due to better cache associativity with most memory accesses to values that are already in the cache(s). Once we have divided the work, Chapel offers lots of means to implement the parallelism but to be a true Sieve of Eratosthenes, we need to have the ability to output the results in order; with many of the convenience mechanisms not doing that, the best/simplest option is likely task parallelism with the output results assigned to an rotary indexed array containing the `sync` results. It turns out that, although the Chapel compiler can sometimes optimize the code so the overhead of creating tasks is not onerous, for this case where the actual tasks are somewhat complex, the compiler can't recognize that an automatically generated thread pool(s) are required so we need to generate the thread pool(s) manually. The code that implements the multi-threading of page segments using thread pools is as follows:

{{works with|Chapel|1.24.1|- compile with the --fast compiler command line flag for full optimization}}

 

type Prime = uint(64);

 

write("Found ", count, " primes up to ", LIMIT);

{{out}}

<pre>The first 25 primes are: 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

 

=={{header|Clojure}}==

(remove (set (mapcat #(range (* % %) n %)

(range 2 (Math/sqrt n))))

 

The above is **not strictly a Sieve of Eratosthenes** as the composite culling ranges (in the ''mapcat'') include all of the multiples of all of the numbers and not just the multiples of primes. When tested with <code clojure>(println (time (count (primes< 1000000))))</code>, it takes about 5.5 seconds just to find the number of primes up to a million, partly because of the extra work due to the use of the non-primes, and partly because of the constant enumeration using sequences with multiple levels of function calls. Although very short, this code is likely only useful up to about this range of a million.

It may be written using the ''into #{}'' function to run slightly faster due to the ''set'' function being concerned with only distinct elements whereas the ''into #{}'' only does the conjunction, and even at that doesn't do that much as it does the conjunction to an empty sequence, the code as follows:

 

(remove (into #{}

(mapcat #(range (* % %) n %)

(range 2 (Math/sqrt n))))

 

The above code is slightly faster for the reasons given, but is still not strictly a Sieve of Eratosthenes due to sieving by all numbers and not just by the base primes.

 

The following code also uses the ''into #{}'' transducer but has been slightly wheel-factorized to sieve odds-only:

(if (< n 2) ()

(cons 2 (remove (into #{}

(mapcat #(range (* % %) n %)

(range 3 (Math/sqrt n) 2)))

 

The above code is a little over twice as fast as the non-odds-only due to the reduced number of operations. It still isn't strictly a Sieve of Eratosthenes as it sieves by all odd base numbers and not only by the base primes.

 

The following code calculates primes up to and including ''n'' using a mutable boolean array but otherwise entirely functional code; it is tens (to a hundred) times faster than the purely functional codes due to the use of mutability in the boolean array:

"Computes lazy sequence of prime numbers up to a given number using sieve of Eratosthenes"

[n]

(do (dorun (map #(cullp %) (filter #(not (aget cmpsts %))

(range 2 (inc root)))))

 

'''Alternative implementation using Clojure's side-effect oriented list comprehension.'''

 

"Returns a lazy sequence of prime numbers less than lim"

[lim]

(doseq [j (range (* i i) lim i)]

(aset refs j false)))

 

'''Alternative implementation using Clojure's side-effect oriented list comprehension. Odds only.'''

"Returns a lazy sequence of prime numbers less than lim"

[lim]

(doseq [j (range (* (+ i i) (inc i)) max-i (+ i i 1))]

(aset refs j false)))

 

This implemantation is about twice as fast as the previous one and uses only half the memory.

'''Alternative very slow entirely functional implementation using lazy sequences'''

 

"Computes lazy sequence of prime numbers up to a given number using sieve of Eratosthenes"

[n]

(cons p (lazy-seq (nxtprm (-> (range (* p p) (inc n) p)

set (remove cs) rest)))))))]

 

The reason that the above code is so slow is that it has has a high constant factor overhead due to using a (hash) set to remove the composites from the future composites stream, each prime composite stream removal requires a scan across all remaining composites (compared to using an array or vector where only the culled values are referenced, and due to the slowness of Clojure sequence operations as compared to iterator/sequence operations in other languages.

 

Here is an immutable boolean vector based non-lazy sequence version other than for the lazy sequence operations to output the result:

"Computes lazy sequence of prime numbers up to a given number using sieve of Eratosthenes"

[max-prime]

(assoc 1 false)

(sieve 2))

 

The above code is still quite slow due to the cost of the immutable copy-on-modify operations.

 

The following code implements an odds-only sieve using a mutable bit packed long array, only using a lazy sequence for the output of the resulting primes:

 

(defn primes-to

(recur (inc i)))))))))]

(if (< n 2) nil

 

The above code is about as fast as any "one large sieving array" type of program in any computer language with this level of wheel factorization other than the lazy sequence operations are quite slow: it takes about ten times as long to enumerate the results as it does to do the actual sieving work of culling the composites from the sieving buffer array. The slowness of sequence operations is due to nested function calls, but primarily due to the way Clojure implements closures by "boxing" all arguments (and perhaps return values) as objects in the heap space, which then need to be "un-boxed" as primitives as necessary for integer operations. Some of the facilities provided by lazy sequences are not needed for this algorithm, such as the automatic memoization which means that each element of the sequence is calculated only once; it is not necessary for the sequence values to be retraced for this algorithm.

 

The following code overcomes many of those limitations by using an internal (OPSeq) "deftype" which implements the ISeq interface as well as the Counted interface to provide immediate count returns (based on a pre-computed total), as well as the IReduce interface which can greatly speed come computations based on the primes sequence (eased greatly using facilities provided by Clojure 1.7.0 and up):

"Computes lazy sequence of prime numbers up to a given number using sieve of Eratosthenes"

[n]

(toString [this] (if (= cnt tcnt) "()"

(.toString (seq (map identity this))))))

 

'(time (count (primes-tox 10000000)))' takes about 40 milliseconds (compiled) to produce 664579.

 

'''A Clojure version of Richard Bird's Sieve using Lazy Sequences (sieves odds only)'''

"Computes the unbounded sequence of primes using a Sieve of Eratosthenes algorithm by Richard Bird."

[]

(do (def oddprms (cons 3 (lazy-seq (let [cmpsts (-> oddprms (allmtpls) (mrgmltpls))]

(minusStrtAt 5 cmpsts)))))

 

The above code is quite slow due to both that the data structure is a linear merging of prime multiples and due to the slowness of the Clojure sequence operations.

 

The following code speeds up the above code by merging the linear sequence of sequences as above by pairs into a right-leaning tree structure:

"Computes the unbounded sequence of primes using a Sieve of Eratosthenes algorithm modified from Bird."

[]

(do (def oddprms (cons 3 (lazy-seq (let [cmpsts (-> oddprms (allmtpls) (mrgmltpls))]

(minusStrtAt 5 cmpsts)))))

 

The above code is still slower than it should be due to the slowness of Clojure's sequence operations.

 

The following code uses a custom "deftype" non-memoizing Co Inductive Stream/Sequence (CIS) implementing the ISeq interface to make the sequence operations more efficient and is about four times faster than the above code:

clojure.lang.ISeq

(first [_] v)

(do (def oddprms (->CIS 3 (fn [] (let [cmpsts (-> oddprms (allmtpls) (mrgmltpls))]

(minusStrtAt 5 cmpsts)))))

 

'(time (count (take-while #(<= (long %) 10000000) (primes-treeFoldingx))))' takes about 3.4 seconds for a range of 10 million.

The following code is a version of the O'Neill Haskell code but does not use wheel factorization other than for sieving odds only (although it could be easily added) and uses a Hash Map (constant amortized access time) rather than a Priority Queue (log n access time for combined remove-and-insert-anew operations, which are the majority used for this algorithm) with a lazy sequence for output of the resulting primes; the code has the added feature that it uses a secondary base primes sequence generator and only adds prime culling sequences to the composites map when they are necessary, thus saving time and limiting storage to only that required for the map entries for primes up to the square root of the currently sieved number:

 

"Infinite sequence of primes using an incremental Sieve or Eratosthenes with a Hashmap"

[]

(cons c (lazy-seq (nxtoddprm (+ c 2) q bsprms cmpsts))))))]

(do (def baseoddprms (cons 3 (lazy-seq (nxtoddprm 5 9 baseoddprms {}))))

 

The above code is slower than the best tree folding version due to the added constant factor overhead of computing the hash functions for every hash map operation even though it has computational complexity of (n log log n) rather than the worse (n log n log log n) for the previous incremental tree folding sieve. It is still about 100 times slower than the sieve based on the bit-packed mutable array due to these constant factor hashing overheads.

In order to implement the O'Neill Priority Queue incremental Sieve of Eratosthenes algorithm, one requires an efficient implementation of a Priority Queue, which is not part of standard Clojure. For this purpose, the most suitable Priority Queue is a binary tree heap based MinHeap algorithm. The following code implements a purely functional (using entirely immutable state) MinHeap Priority Queue providing the required functions of (emtpy-pq) initialization, (getMin-pq pq) to examinte the minimum key/value pair in the queue, (insert-pq pq k v) to add entries to the queue, and (replaceMinAs-pq pq k v) to replaace the minimum entry with a key/value pair as given (it is more efficient that if functions were provided to delete and then re-insert entries in the queue; there is therefore no "delete" or other queue functions supplied as the algorithm does not requrie them:

 

Object

(toString [_] (str "<" k "," v ">")))

(if (<= kl kr)

(recur l #(cont (->PQNode kvl % r)))

 

Note that the above code is written partially using continuation passing style so as to leave the "recur" calls in tail call position as required for efficient looping in Clojure; for practical sieving ranges, the algorithm could likely use just raw function recursion as recursion depth is unlikely to be used beyond a depth of about ten or so, but raw recursion is said to be less code efficient.

The actual incremental sieve using the Priority Queue is as follows, which code uses the same optimizations of postponing the addition of prime composite streams to the queue until the square root of the currently sieved number is reached and using a secondary base primes stream to generate the primes composite stream markers in the queue as was used for the Hash Map version:

 

"Infinite sequence of primes using an incremental Sieve or Eratosthenes with a Priority Queue"

[]

(cons c (lazy-seq (nxtoddprm (+ c 2) q bsprms cmpsts)))))))]

(do (def baseoddprms (cons 3 (lazy-seq (nxtoddprm 5 9 baseoddprms (empty-pq)))))

 

The above code is faster than the Hash Map version up to about a sieving range of fifteen million or so, but gets progressively slower for larger ranges due to having (n log n log log n) computational complexity rather than the (n log log n) for the Hash Map version, which has a higher constant factor overhead that is overtaken by the extra "log n" factor.

To show that Clojure does not need to be particularly slow, the following version runs about twice as fast as the non-segmented unbounded array based version above (extremely fast compared to the non-array based versions) and only a little slower than other equivalent versions running on virtual machines: C# or F# on DotNet or Java and Scala on the JVM:

 

 

(def PGSZ (bit-shift-left 1 14)) ;; size of CPU cache

(next pgseq)

(+ cnt (count-pg PGBTS pg))))))]

 

The above code runs just as fast as other virtual machine languages when run on a 64-bit JVM; however, when run on a 32-bit JVM it runs almost five times slower. This is likely due to Clojure only using 64-bit integers for integer operations and these operations getting JIT compiled to use library functions to simulate those operations using combined 32-bit operations under a 32-bit JVM whereas direct CPU operations can be used on a 64-bit JVM

 

The base primes culling page size is reduced from the page size for the main primes so that there is less overhead for smaller primes ranges; otherwise excess base primes are generated for fairly small sieve ranges.

 

=={{header|CMake}}==

# Check for integer overflow. With CMake using 32-bit signed integer,

# this check fails when limit > 46340.

endforeach(i)

set(${var} ${list} PARENT_SCOPE)

# Print all prime numbers through 100.

eratosthenes(primes 100)

 

=={{header|COBOL}}==

*> which are kludges to get syntax highlighting to work.

IDENTIFICATION DIVISION.

 

GOBACK

 

=={{header|Comal}}==

{{trans|BASIC}}

input "Limit? ": limit

dim sieve(1:limit)

endif

endfor i

 

{{Out}}

 

=={{header|Common Lisp}}==

(loop

with sieve = (make-array (1+ maximum)

and do (loop for composite from (expt candidate 2)

to maximum by candidate

 

Working with odds only (above twice speedup), and marking composites only for primes up to the square root of the maximum:

 

"Prime numbers sieve for odd numbers.

Returns a list with all the primes that are less than or equal to maximum."

:do (loop :for j :from (* i (1+ i) 2) :to maxi :by odd-number

:do (setf (sbit sieve j) 1))

 

The indexation scheme used here interprets each index <code>i</code> as standing for the value <code>2i+1</code>. Bit <code>0</code> is unused, a small price to pay for the simpler index calculations compared with the <code>2i+3</code> indexation scheme. The multiples of a given odd prime <code>p</code> are enumerated in increments of <code>2p</code>, which corresponds to the index increment of <code>p</code> on the sieve array. The starting point <code>p*p = (2i+1)(2i+1) = 4i(i+1)+1</code> corresponds to the index <code>2i(i+1)</code>.

 

=={{header|Cowgol}}==

 

# To change the maximum prime, change the size of this array

end if;

cand := cand + 1;

 

{{out}}

4987

4999</pre>

 

=={{header|Crystal}}==

This implementation uses a `BitArray` so it is automatically bit-packed to use just one bit per number representation:

 

 

require "bit_array"

print SoE.new(1_000_000_000).each.size

elpsd = (Time.monotonic - start_time).total_milliseconds

 

{{out}}

the non-odds-only version as per the above should never be used because in not using odds-only, it uses twice the memory and over two and a half times the CPU operations as the following odds-only code, which is very little more complex:

 

 

require "bit_array"

print SoE_Odds.new(1_000_000_000).each.size

elpsd = (Time.monotonic - start_time).total_milliseconds

 

{{out}}

For sieving of ranges larger than a few million efficiently, a page-segmented sieve should always be used to preserve CPU cache associativity by making the page size to be about that of the CPU L1 data cache. The following code implements a page-segmented version that is an extensible sieve (no upper limit needs be specified) using a secondary memoized feed of base prime value arrays which use a smaller page-segment size for efficiency. When the count of the number of primes is desired, the sieve is polymorphic in output and counts the unmarked composite bits by using fast `popcount` instructions taken 64-bits at a time. The code is as follows:

 

 

alias Prime = UInt64

answr = primes_count_to(1_000_000_000) # fast way

elpsd = (Time.monotonic - start_time).total_milliseconds

 

{{out}}

=={{header|D}}==

===Simpler Version===

 

uint[] sieve(in uint limit) nothrow @safe {

void main() {

50.sieve.writeln;

{{out}}

<pre>[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]</pre>

===Faster Version===

This version uses an array of bits (instead of booleans, that are represented with one byte), and skips even numbers. The output is the same.

 

size_t[] sieve(in size_t m) pure nothrow @safe {

void main() {

50.sieve.writeln;

 

===Extensible Version===

(This version is used in the task [[Extensible prime generator#D|Extensible prime generator]].)

struct Prime {

uint[] a = [2];

uint.max.iota.map!prime.until!q{a > 50}.writeln;

}

To see the output (that is the same), compile with <code>-version=sieve_of_eratosthenes3_main</code>.

 

=={{header|Dart}}==

String iterableToString(Iterable seq) {

String str = "[";

print(iterableToString(sortedValues)); // expect sieve.length to be 168 up to 1000...

// Expect.equals(168, sieve.length);

{{out}}<pre>

Found 168 primes up to 1000 in 9 milliseconds.

 

===faster bit-packed array odds-only solution===

import 'dart:math';

 

print("There were $answer primes found up to $range.");

print("This test bench took $elapsed milliseconds.");

{{output}}

<pre>( 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 )

The following code will have about O(n log (log n)) performance due to a hash table having O(1) average performance and is only somewhat slow due to the constant overhead of processing hashes:

 

Iterable<int> oddprms() sync* {

yield(3); yield(5); // need at least 2 for initialization

print("There were $answer primes found up to $range.");

print("This test bench took $elapsed milliseconds.");

{{output}}

<pre>( 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 )

 

{{trans|Kotlin}}

import 'dart:math';

import 'dart:collection';

print("There were $answer primes found up to $range.");

print("This test bench took $elapsed milliseconds.");

{{output}}

<pre>( 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 )

 

The algorithm can be sped up by a factor of four by extreme wheel factorization and (likely) about a factor of the effective number of CPU cores by using multi-processing isolates, but there isn't much point if one is to use the prime generator for output. For most purposes, it is better to use custom functions that directly manipulate the culled bit-packed page segments as `countPrimesTo` does here.

 

=={{header|Delphi}}==

 

{$APPTYPE CONSOLE}

Sieve.Free;

readln;

Output:

<pre>2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 </pre>

 

=={{header|DWScript}}==

 

var

n, k : Integer;

var i : Integer;

for i:=0 to r.High do

 

=={{header|Dylan}}==

With outer to sqrt and inner to p^2 optimizations:

let limit = floor(n ^ 0.5) + 1;

let sieve = make(limited(<simple-vector>, of: <boolean>), size: n + 1, fill: #t);

if (sieve[x]) format-out("Prime: %d\n", x); end;

end;

 

=={{header|E}}==

 

=={{header|EasyLang}}==

.

 

=={{header|eC}}==

{{incorrect|eC|It uses rem testing and so is a trial division algorithm, not a sieve of Eratosthenes.}}

Note: this is not a Sieve of Eratosthenes; it is just trial division.

public class FindPrime

{

}

}

 

=={{header|EchoLisp}}==

===Sieve===

 

;; converts sieve->list for integers in [nmin .. nmax[

→ (1000003 1000033 1000037 1000039 1000081 1000099)

(s-next-prime s-primes 9_000_000)

 

===Segmented sieve===

Allow to extend the basis sieve (n) up to n^2. Memory requirement is O(√n)

;; delta multiple of sqrt(n)

;; segment is [left .. left+delta-1]

 

;; 8 msec using the native (prime?) function

 

===Wheel===

A 2x3 wheel gives a 50% performance gain.

(define (weratosthenes n)

(define primes (make-bit-vector n )) ;; everybody to #f (false)

(for ([j (in-range (* p p) n p)])

(bit-vector-set! primes j #f)))

 

=={{header|EDSAC order code}}==

 

On the EdsacPC simulator (see link above) the printout starts off very slowly, and gradually gets faster.

[Sieve of Eratosthenes]

[EDSAC program. Initial Orders 2]

E 19 Z [address to start execution]

P F [acc = 0 at start]

{{out}}

<pre>

=={{header|Eiffel}}==

{{works with|EiffelStudio|6.6 beta (with provisional loop syntax)}}

APPLICATION

end

end

 

Output:

 

=={{header|Elixir}}==

def eratosthenes(limit \\ 1000) do

sieve = [false, false | Enum.to_list(2..limit)] |> List.to_tuple

if x=elem(sieve, n), do: :io.format("~3w", [x]), else: :io.format(" .")

if rem(n+1, 20)==0, do: IO.puts ""

 

{{out}}

Shorter version (but slow):

 

defmodule Sieve do

def primes_to(limit), do: sieve(Enum.to_list(2..limit))

defp sieve([]), do: []

end

 

'''Alternate much faster odds-only version more suitable for immutable data structures using a (hash) Map'''

 

The above code has a very limited useful range due to being very slow: for example, to sieve to a million, even changing the algorithm to odds-only, requires over 800 thousand "copy-on-update" operations of the entire saved immutable tuple ("array") of 500 thousand bytes in size, making it very much a "toy" application. The following code overcomes that problem by using a (immutable/hashed) Map to store the record of the current state of the composite number chains resulting from each of the secondary streams of base primes, which are only 167 in number up to this range; it is a functional "incremental" Sieve of Eratosthenes implementation:

@typep stt :: {integer, integer, integer, Enumerable.integer, %{integer => integer}}

 

|> (fn {t,ans} ->

IO.puts "There are #{ans} primes up to #{range}."

{{output}}

<pre>The first 25 primes are:

 

In order to save the computation time of computing the hashes, the following version uses a deferred execution Co-Inductive Stream type (constructed using Tuple's) in an infinite tree folding structure (by the `pairs` function):

@typep cis :: {integer, (() -> cis)}

@typep ciss :: {cis, (() -> ciss)}

|> (fn {t,ans} ->

IO.puts "There are #{ans} primes up to #{range}."

 

It's output is identical to the previous version other than the time required is less than half; however, it has a O(n (log n) (log (log n))) asymptotic computation complexity meaning that it gets slower with range faster than the above version. That said, it would take sieving to billions taking hours before the two would take about the same time.

 

=={{header|Emacs Lisp}}==

(let ((xs (make-vector (1+ limit) 0)))

 

Straightforward implementation of [http://en.wikipedia.org/wiki/Sieve_of_Eratosthenes#Implementation sieve of Eratosthenes], 2 times faster:

 

(let ((xs (vconcat [0 0] (number-sequence 2 limit))))

 

=={{header|Erlang}}==

===Erlang using Dicts===

-module( sieve_of_eratosthenes ).

 

 

find_prime( error, _N, Acc ) -> Acc;

{{out}}

<pre>

Has the virtue of working for any -> N :)

 

-module( sieve ).

-export( [main/1,primes/2] ).

Primes = lists:filter( fun({_,Y}) -> Y > 0 end, Tuples),

[ X || {X,_} <- Primes ].

{{out}}

<pre>

Since I had written a really odd and slow one, I thought I'd best do a better performer. Inspired by an example from https://github.com/jupp0r

 

 

-module(ossieve).

ResultSet = ordsets:add_element(2,sieve(Candidates,Candidates,ordsets:new(),N)),

io:fwrite("Sieved... ~w~n",[ResultSet]).

{{out}}

<pre>

A pure list comprehension approach.

 

-module(sieveof).

-export([main/1,primes/1, primes/2]).

remove([H | X], [H | Y]) -> remove(X, Y);

remove(X, [H | Y]) -> [H | remove(X, Y)].

{out}

<pre>

===Erlang ets + cpu distributed implementation ===

much faster previous erlang examples

#!/usr/bin/env escript

%% -*- erlang -*-

comp_i(J, I, N) when J =< N -> ets:insert(comp, {J, 1}), comp_i(J+I, I, N);

comp_i(J, _, N) when J > N -> ok.

{{out}}

<pre>

 

another several erlang implementation: http://mijkenator.github.io/2015/11/29/project-euler-problem-10/

 

=={{header|ERRE}}==

PROGRAM SIEVE_ORG

! --------------------------------------------------

PRINT(COUNT%;" PRIMES")

END PROGRAM

{{out}}

last lines of the output screen

16301 16319 16333 16339 16349 16361 16363 16369 16381

1899 PRIMES

</pre>

 

=={{header|Euphoria}}==

sequence flags,primes

flags = repeat(1, limit)

end if

end for

 

Output:

=={{header|F Sharp}}==

===Short with mutable state===

let primes max =

let mutable xs = [|2..max|]

xs <- xs |> Array.except [|x*x..x..max|]

xs

===Short Sweet Functional and Idiotmatic===

Well lists may not be lazy, but if you call it a sequence then it's a lazy list!

(*

An interesting implementation of The Sieve of Eratosthenes.

yield g; yield! fn (g - 1) g |> Seq.map2 (&&) ng |> Seq.cache |> fg }

Seq.initInfinite (fun x -> true) |> fg

{{out}}

<pre>

 

This is the idea behind Richard Bird's unbounded code presented in the Epilogue of [http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf M. O'Neill's article] in Haskell. It is about twice as much code as the Haskell code because F# does not have a built-in lazy list so that the effect must be constructed using a Co-Inductive Stream (CIS) type since no memoization is required, along with the use of recursive functions in combination with sequences. The type inference needs some help with the new CIS type (including selecting the generic type for speed). Note the use of recursive functions to implement multiple non-sharing delayed generating base primes streams, which along with these being non-memoizing means that the entire primes stream is not held in memory as for the original Bird code:

 

let primesBird() =

else minusat (n + 1u) (ctlf())

let rec baseprms() = CIS(2u, fun() -> baseprms() |> allmltps |> cmpsts |> minusat 3u)

 

The above code sieves all numbers of two and up including all even numbers as per the page specification; the following code makes the very minor changes for an odds-only sieve, with a speedup of over a factor of two:

 

let primesBirdOdds() =

else minusat (n + 2u) (ctlf())

let rec oddprms() = CIS(3u, fun() -> oddprms() |> allmltps |> cmpsts |> minusat 5u)

 

'''Tree Folding Sieve'''

 

The above code is still somewhat inefficient as it operates on a linear right extending structure that deepens linearly with increasing base primes (those up to the square root of the currently sieved number); the following code changes the structure into an infinite binary tree-like folding by combining each pair of prime composite streams before further processing as usual - this decreases the processing by approximately a factor of log n:

 

let primesTreeFold() =

else minusat (n + 2u) (ctlf())

let rec oddprms() = CIS(3u, fun() -> oddprms() |> allmltps |> cmpsts |> minusat 5u)

 

The above code is over four times faster than the "BirdOdds" version (at least 10x faster than the first, "primesBird", producing the millionth prime) and is moderately useful for a range of the first million primes or so.

 

In order to investigate Priority Queue Sieves as espoused by O'Neill in the referenced article, one must find an equivalent implementation of a Min Heap Priority Queue as used by her. There is such an purely functional implementation [http://rosettacode.org/wiki/Priority_queue#Functional in RosettaCode translated from the Haskell code she used], from which the essential parts are duplicated here (Note that the key value is given an integer type in order to avoid the inefficiency of F# in generic comparison):

module MinHeap =

 

 

let replaceMin wk wv = function | Mt -> Mt

 

Except as noted for any individual code, all of the following codes need the following prefix code in order to implement the non-memoizing Co-Inductive Streams (CIS's) and to set the type of particular constants used in the codes to the same time as the "Prime" type:

type Prime = uint32

let frstprm = 2u

let frstoddprm = 3u

let inc1 = 1u

 

The F# equivalent to O'Neill's "odds-only" code is then implemented as follows, which needs the included changed prefix in order to change the primes type to a larger one to prevent overflow (as well the key type for the MinHeap needs to be changed from uint32 to uint64); it is functionally the same as the O'Neill code other than for minor changes to suit the use of CIS streams and the option output of the "peekMin" function:

type Prime = uint64

let frstprm = 2UL

let rec nxto i = CIS(i, fun() -> nxto (i + inc)) in nxto frstoddprm

Seq.unfold (fun (cis: CIS<Prime>) -> Some(cis.v, cis.cont()))

 

However, that algorithm suffers in speed and memory use due to over-eager adding of prime composite streams to the queue such that the queue used is much larger than it needs to be and a much larger range of primes number must be used in order to avoid numeric overflow on the square of the prime added to the queue. The following code corrects that by using a secondary (actually a multiple of) base primes streams which are constrained to be based on a prime that is no larger than the square root of the currently sieved number - this permits the use of much smaller Prime types as per the default prefix:

let rec nxtprm n pq q (bps: CIS<Prime>) =

if n >= q then let bp = bps.v in let adv = bp + bp

nxtprm (frstoddprm + inc) MinHeap.empty (frstoddprm * frstoddprm) (oddprms()))

Seq.unfold (fun (cis: CIS<Prime>) -> Some(cis.v, cis.cont()))

 

The above code is well over five times faster than the previous translated O'Neill version for the given variety of reasons.

The following code is written in functional style other than it uses a mutable bit array to sieve the composites:

 

let buf = System.Collections.BitArray(int limit + 1, true)

let cull p = { p * p .. p .. limit } |> Seq.iter (fun c -> buf.[int c] <- false)

if argv = null || argv.Length = 0 then failwith "no command line argument for limit!!!"

printfn "%A" (primes (System.UInt32.Parse argv.[0]) |> Seq.length)

 

Substituting the following minor changes to the code for the "primes" function will only deal with the odd prime candidates for a speed up of over a factor of two as well as a reduction of the buffer size by a factor of two:

 

let lmtb,lmtbsqrt = (limit - 3u) / 2u, (uint32 (sqrt (double limit)) - 3u) / 2u

let buf = System.Collections.BitArray(int lmtb + 1, true)

let oddprms = { 0u .. lmtb } |> Seq.map (fun i -> if buf.[int i] then i + i + 3u else 0u)

|> Seq.filter ((<>) 0u)

 

The following code uses other functional forms for the inner culling loops of the "primes function" to reduce the use of inefficient sequences so as to reduce the execution time by another factor of almost three:

 

let lmtb,lmtbsqrt = (limit - 3u) / 2u, (uint32 (sqrt (double limit)) - 3u) / 2u

let buf = System.Collections.BitArray(int lmtb + 1, true)

let rec cullp c = if c <= lmtb then buf.[int c] <- false; cullp (c + p)

(if buf.[int i] then cullp s); culltest (i + 1u) in culltest 0u

 

Now much of the remaining execution time is just the time to enumerate the primes as can be seen by turning "primes" into a primes counting function by substituting the following for the last line in the above code doing the enumeration; this makes the code run about a further five times faster:

 

if i > int lmtb then acc else if buf.[i] then count (i + 1) (acc + 1) else count (i + 1) acc

 

Since the final enumeration of primes is the main remaining bottleneck, it is worth using a "roll-your-own" enumeration implemented as an object expression so as to save many inefficiencies in the use of the built-in seq computational expression by substituting the following code for the last line of the previous codes, which will decrease the execution time by a factor of over three (instead of almost five for the counting-only version, making it almost as fast):

 

let i = ref -2

let rec nxti() = i:=!i + 1;if !i <= int lmtb && not buf.[!i] then nxti() else !i <= int lmtb

member this.GetEnumerator() = nmrtr()

interface System.Collections.IEnumerable with

 

The various optimization techniques shown here can be used "jointly and severally" on any of the basic versions for various trade-offs between code complexity and performance. Not shown here are other techniques of making the sieve faster, including extending wheel factorization to much larger wheels such as 2/3/5/7, pre-culling the arrays, page segmentation, and multi-processing.

 

The following code uses the DotNet Dictionary class instead of the above functional Priority Queue to implement the sieve; as average (amortized) hash table access is O(1) rather than O(log n) as for the priority queue, this implementation is slightly faster than the priority queue version for the first million primes and will always be faster for any range above some low range value:

let frstprm = 2u

let frstoddprm = 3u

nxtprm (frstoddprm + inc) (frstoddprm * frstoddprm) (oddprms()))

Seq.unfold (fun (cis: CIS<Prime>) -> Some(cis.v, cis.cont()))

 

The above code uses functional forms of code (with the imperative style commented out to show how it could be done imperatively) and also uses a recursive non-sharing secondary source of base primes just as for the Priority Queue version. As for the functional codes, the Primes type can easily be changed to "uint64" for wider range of sieving.

 

All of the above unbounded implementations including the above Dictionary based version are quite slow due to their large constant factor computational overheads, making them more of an intellectual exercise than something practical, especially when larger sieving ranges are required. The following code implements an unbounded page segmented version of the sieve in not that many more lines of code, yet runs about 25 times faster than the Dictionary version for larger ranges of sieving such as to one billion; it uses functional forms without mutability other than for the contents of the arrays and the `primes` enumeration generator function that must use mutability for speed:

type private PrimeNdx = float // they are slow in JavaScript polyfills

 

let answr = countPrimesTo limit // over twice as fast way!

let elpsd = (System.DateTime.Now.Ticks - start) / 10000L

 

{{out}}

 

Factor is pleasantly multiparadigm. Usually, it's natural to write more functional or declarative code in Factor, but this is an instance where it is more natural to write imperative code. Lexical variables are useful here for expressing the necessary mutations in a clean way.

math.ranges prettyprint sequences ;

IN: rosetta-code.sieve-of-erato

120 sieve . ;

 

 

=={{header|FOCAL}}==

1.2 A N

1.3 I (2047-N)5.1

4.1 T %4.0,W,!

 

Note that with the 4k paper tape version of FOCAL, the program will run out of memory for N>190 or so.

 

The above code is not really very good Forth style as the main initialization, sieving, and output, are all in one `sieve` routine which makes it difficult to understand and refactor; Forth code is normally written in a series of very small routines which makes it easier to understand what is happening on the data stack, since Forth does not have named local re-entrant variable names as most other languages do for local variables (which other languages also normally store local variables on the stack). Also, it uses the `HERE` pointer to user space which points to the next available memory after all compilation is done as a unsized buffer pointer, but as it does not reserve that space for the sieving buffer, it can be changed by other concatenated routines in unexpected ways; better is to allocate the sieving buffer as required from the available space at the time the routines are run and pass that address between concatenated functions until a finalization function frees the memory and clears the stack; this is equivalent to allocating from the "heap" in other languages. The below code demonstrates these ideas:

 

 

\ given square index and prime index, u0, sieve the multiples of said prime...

\ testing the code...

100 initsieve sieve .primes

 

{{out}}

Although the above version resolves many problems of the first version, it is wasteful of memory as each composite number in the sieve buffer is a byte of eight bits representing a boolean value. The memory required can be reduced eight-fold by bit packing the sieve buffer; this will take more "bit-twiddling" to read and write the bits, but reducing the memory used will give better cache assiciativity to larger ranges such that there will be a net gain in performance. This will make the code more complex and the stack manipulations will be harder to write, debug, and maintain, so ANS Forth 1994 provides a local variable naming facility to make this much easier. The following code implements bit-packing of the sieve buffer using local named variables when required:

 

: numbts ( u -- u ) \ pop count number of bits...

0 SWAP BEGIN DUP 0<> WHILE SWAP 1+ SWAP DUP 1- AND REPEAT DROP ;

 

100 initsieve sieve .primes

 

The output of the above code is the same as the previous version, but it takes about two thirds the time while using eight times less memory; it takes about 6.5 seconds on my Intel Skylake i5-6500 at 3.6 GHz (turbo) using swiftForth (32-bit) and about 3.5 seconds on VFX Forth (64-bit), both of which compile to machine code but with the latter much more optimized; gforth-fast is about twice as slow as swiftForth and five times slower then VFX Forth as it just compiles to threaded execution tokens (more like an interpreter).

While the above version does greatly reduce the amount of memory used for a given sieving range and thereby also somewhat reduces execution time; any sieve intended for sieving to limits of a hundred million or more should use a page-segmented implementation; page-segmentation means that only storage for a representation of the base primes up to the square root of the limit plus a sieve buffer that should also be at least proportional to the same square root is required; this will again make the execution faster as ranges go up due to better cache associativity with most memory accesses being within the CPU cache sizes. The following Forth code implements a basic version that does this:

 

1 17 LSHIFT CONSTANT L1CacheBits

L1CacheBits 8 * CONSTANT L2CacheBits

 

100 .primes

 

{{out}}

 

=={{header|Fortran}}==

{{works with|Fortran|90 and later}}

 

implicit none

write (*, *)

 

Output:

 

Because it uses four byte logical's (default size) as elements of the sieve buffer, the above code uses 400 bytes of memory for this trivial task of sieving to 100; it also has 49 + 31 + 16 + 8 = 104 (for the culling by the primes of two, three, five, and seven) culling operations.

'''Optimised using a pre-computed wheel based on 2:'''

 

 

implicit none

write (*, *)

 

Output:

 

This so-called "optimized" version still uses 400 bytes of memory but slightly reduces to 74 operations from 104 operations including the initialization of marking all of the even representations as composite due to skipping the re-culling of the even representation, so isn't really much of an optimization at all!

The above implementations, especially the second odds-only code, are some of the most inefficient versions of the Sieve of Eratosthenes in any language here as to time and space efficiency, only worse by some naive JavaScript implementations that use eight-byte Number's as logical values; the second claims to be wheel factorized but still uses all the same memory as the first and still culls by the even numbers in the initialization of the sieve buffer. As well, using four bytes (default logical size) to store a boolean value is terribly wasteful if these implementations were to be extended to non-toy ranges. The following code implements proper wheel factorization by two, reducing the space used by a factor of about eight to 49 bytes by using `byte` as the sieve buffer array elements and not requiring the evens initialization, thus reducing the number of operations to 16 + 8 + 4 = 28 (for the culling primes of three, five, and seven) culling operations:

 

implicit none

write (*, *)

 

The output is the same as the earlier version.

The above implementation is still space inefficient in effectively only using one bit out of eight; the following version implements bit packing to reduce memory use by a factor of eight by using bits to represent composite numbers rather than bytes:

 

implicit none

print '(a, i0, a, i0, a, f0.0, a)', &

'There are ', cnt, ' primes up to ', i_max, '.'

 

{{out}}

As well as adding page-segmentation, the following code adds multi-processing which is onc of the capabilities for which modern Fortran is known:

 

 

implicit none

deallocate(sbaa)

 

{{out}}

===Basic version===

function Sieve returns a list of primes less than or equal to the given aLimit

program prime_sieve;

{$mode objfpc}{$coperators on}

I: DWord;

begin

aList.Free;

end;

end;

end.

===Alternative segmented(odds only) version===

function OddSegmentSieve returns a list of primes less than or equal to the given aLimit

program prime_sieve;

{$mode objfpc}{$coperators on}

for I := 0 to Length(aList) - 2 do

Write(aList[I], ', ');

end;

begin

PrintPrimes(OddSegmentSieve(Limit));

end.

 

=={{header|FreeBASIC}}==

 

Sub sieve(n As Integer)

Print

Print "Press any key to quit"

 

{{out}}

 

=={{header|Frink}}==

n = eval[input["Enter highest number: "]]

results = array[sieve[n]]

return select[array, { |x| x != 0 }]

}

 

=={{header|Furor}}==

''Note: With benchmark function''

 

tick sto startingtick

#g 100000 sto MAX

end

{ „MAX” } { „startingtick” } { „primeNumbers” } { „count” }

 

 

 

=={{header|FutureBasic}}==

===Basic sieve of array of booleans===

 

begin globals

end globals

 

local fn SieveOfEratosthenes( n as long )

end fn

 

fn SieveOfEratosthenes( 100 )

Output:

<pre>

=={{header|Fōrmulæ}}==

 

 

 

 

=={{header|GAP}}==

local a, i, j;

a := ListWithIdenticalEntries(n, true);

Eratosthenes(100);

 

 

=={{header|GLBasic}}==

// GLBasic implementation

 

 

KEYWAIT

 

=={{header|Go}}==

===Basic sieve of array of booleans===

import "fmt"

 

}

}

Output:

<pre>

The above version's output is rather specialized; the following version uses a closure function to enumerate over the culled composite number array, which is bit packed. By using this scheme for output, no extra memory is required above that required for the culling array:

 

 

import (

}

fmt.Printf("\r\n%v\r\n", count)

{{output}}

<pre>2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

===Sieve Tree===

A fairly odd sieve tree method:

import "fmt"

 

fmt.Println(p())

}

 

===Concurrent Daisy-chain sieve===

A concurrent prime sieve adopted from the example in the "Go Playground" window at http://golang.org/

import "fmt"

}

}

The output:

<pre>

===Postponed Concurrent Daisy-chain sieve===

Here we postpone the ''creation'' of filters until the prime's square is seen in the input, to radically reduce the amount of filter channels in the sieve chain.

import "fmt"

}

}

 

The output:

===Incremental Odds-only Sieve===

Uses Go's built-in hash tables to store odd composites, and defers adding new known composites until the square is seen.

package main

 

 

var p int

p = <-primes

 

}

}

The output:

<pre>

=={{header|Groovy}}==

This solution uses a BitSet for compactness and speed, but in [[Groovy]], BitSet has full List semantics. It also uses both the "square root of the boundary" shortcut and the "square of the prime" shortcut.

def isPrime = new BitSet(bound)

isPrime[0..1] = false

}

(0..bound).findAll { isPrime[it] }

 

Test:

 

Output:

=={{header|GW-BASIC}}==

 

20 DIM FLAGS(LIMIT)

30 FOR N = 2 TO SQR (LIMIT)

100 FOR N = 2 TO LIMIT

110 IF FLAGS(N) = 0 THEN PRINT N;", ";

 

=={{header|Haskell}}==

Mutable array of unboxed <code>Bool</code>s indexed by <code>Int</code>s:

 

{-# OPTIONS_GHC -O2 #-}

 

 

putStrLn $ "Found " ++ show answr ++ " to " ++ show top ++

 

The above code chooses conciseness and elegance over speed, but it isn't too slow:

Mutable array of unboxed <code>Bool</code>s indexed by <code>Int</code>s, representing odds only:

 

import Control.Monad.ST

import Data.Array.ST

primesToUO :: Int -> [Int]

primesToUO top | top > 1 = 2 : [2*i + 1 | (i,True) <- assocs $ sieveUO top]

 

This represents ''odds only'' in the array. [http://ideone.com/KwZNc Empirical orders of growth] is ~ <i>n^1.2</i> in ''n'' primes produced, and improving for bigger ''n''&zwj;&thinsp;&zwj;s. Memory consumption is low (array seems to be packed) and growing about linearly with ''n''. Can further be [http://ideone.com/j24jxV significantly sped up] by re-writing the <code>forM_</code> loops with direct recursion, and using <code>unsafeRead</code> and <code>unsafeWrite</code> operations.

The reason for this alternate version is to have an accessible version of "odds only" that uses the same optimizations and is written in the same coding style as the basic version. This can be used by just substituting the following code for the function of the same name in the first base example above. Mutable array of unboxed <code>Bool</code>s indexed by <code>Int</code>s, representing odds only:

 

primesTo limit

| limit < 2 = []

$ takeWhile ((>=) lmt . snd) -- for bp's <= square root limit

[ getbpndx i | (i, False) <- assocs oddcmpstsf ]

 

{{out}}

===Immutable arrays===

Monolithic sieving array. ''Even'' numbers above 2 are pre-marked as composite, and sieving is done only by ''odd'' multiples of ''odd'' primes:

primesToA m = sieve 3 (array (3,m) [(i,odd i) | i<-[3..m]] :: UArray Int Bool)

| p*p > m = 2 : [i | (i,True) <- assocs a]

| a!p = sieve (p+2) $ a//[(i,False) | i <- [p*p, p*p+2*p..m]]

 

Its performance sharply depends on compiler optimizations. Compiled with -O2 flag in the presence of the explicit type signature, it is very fast in producing first few million primes. <code>(//)</code> is an array update operator.

 

Works by segments between consecutive primes' squares. Should be the fastest of non-monadic code. ''Evens'' are entirely ignored:

 

primesSA = 2 : prs ()

a :: UArray Int Bool

a = accumArray (\ b c -> False) True (1,q-1)

 

===Basic list-based sieve===

Straightforward implementation of the sieve of Eratosthenes in its original bounded form. This finds primes in gaps between the composites, and composites as an enumeration of each prime's multiples.

eratos (p : xs)

| p*p > m = p : xs

EQ -> minus xs ys

GT -> minus a ys

Its time complexity is similar to that of optimal [[Primality_by_trial_division#Haskell|trial division]] because of limitations of Haskell linked lists, where <code>(minus a b)</code> takes time proportional to <code>length(union a b)</code> and not <code>(length b)</code>, as achieved in imperative setting with direct-access memory. Uses ordered list representation of sets.

 

===Unbounded list based sieve===

Unbounded, "naive", too eager to subtract (see above for the definition of <code>minus</code>):

where

sieve (p:xs) = p : sieve (minus xs [p, p+p..])

This is slow, with complexity increasing as a square law or worse so that it is only moderately useful for the first few thousand primes or so.

 

The number of active streams can be limited to what's strictly necessary by postponement until the square of a prime is seen, getting a massive complexity improvement to better than <i>~ n^1.5</i> so it can get first million primes or so in a tolerable time:

where

sieve (x:xs) q (p:t)

-- fix $ (2:) . concat

-- . unfoldr (\(p:ps,xs)-> Just . second ((ps,) . (`minus` [p*p, p*p+p..]))

 

Transposing the workflow, going by segments between the consecutive squares of primes:

 

primesSE = 2 : sieve 3 4 (tail primesSE) (inits primesSE)

-- True (x,q-1) [(i,()) | f <- fs, let n=div(x+f-1)f*f,

-- i <- [n, n+f..q-1]] :: UArray Int Bool )]

 

The basic gradually-deepening left-leaning <code>(((a-b)-c)- ... )</code> workflow of <code>foldl minus a bs</code> above can be rearranged into the right-leaning <code>(a-(b+(c+ ... )))</code> workflow of <code>minus a (foldr union [] bs)</code>. This is the idea behind Richard Bird's unbounded code presented in [http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf M. O'Neill's article], equivalent to:

 

 

-- = _Y ( (2:) . minus [3..] . _LU . map(\p-> [p*p, p*p+p..]) )

LT -> x : union xs b

EQ -> x : union xs ys

 

Using <code>_Y</code> is meant to guarantee the separate supply of primes to be independently calculated, recursively, instead of the same one being reused, corecursively; thus the memory footprint is drastically reduced. This idea was introduced by M. ONeill as a double-staged production, with a separate primes feed.

 

This merges primes' multiples streams in a ''tree''-like fashion, as a sequence of balanced trees of <code>union</code> nodes, likely achieving theoretical time complexity only a ''log n'' factor above the optimal ''n log n log (log n)'', for ''n'' primes produced. Indeed, empirically it runs at about ''~ n^1.2'' (for producing first few million primes), similarly to priority-queue&ndash;based version of M. O'Neill's, and with very low space complexity too (not counting the produced sequence of course):

primes() = 2 : _Y ((3:) . gaps 5 . _U . map(\p-> [p*p, p*p+2*p..])) where

_Y g = g (_Y g) -- = g (g (g ( ... ))) non-sharing multistage fixpoint combinator

_U ((x:xs):t) = x : (merge xs . _U . pairs) t -- tree-shaped folding big union

pairs (xs:ys:t) = merge xs ys : pairs t

 

Works with odds only, the simplest kind of wheel. Here's the [http://ideone.com/qpnqe test entry] on Ideone.com, and a [http://ideone.com/p0e81 comparison with more versions].

====With Wheel====

Using <code>_U</code> defined above,

primesW = [2,3,5,7] ++ _Y ( (11:) . gapsW 13 (tail wheel) . _U .

map (\p->

wheel = 2:4:2:4:6:2:6:4:2:4:6:6:2:6:4:2:6:4:6:8:4:2:4:2: -- gaps = (`gapsW` cycle [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:wheel

 

Used [[Emirp_primes#List-based|here]] and [[Extensible_prime_generator#List_based|here]].

2. Improving the means to re-generate the position on the wheel for the recursive base primes without the use of `dropWhile`, etc. The below improved code uses a copy of the place in the wheel for each found base prime for ease of use in generating the composite number to-be-culled chains.

 

wheelGen :: Int -> ([Int],Int,[Int])

wheelGen n = loop 1 3 [2] [2] where

| x < y = x : union xs' ys

| y < x = y : union xs ys'

 

When compiled with -O2 optimization and -fllvm (the LLVM back end), the above code is over twice as fast as the Odds-Only version as it should be as that is about the ratio of reduced operations minus some slightly increased operation complexity, sieving the primes to a hundred million in about seven seconds on a modern middle range desktop computer. It is almost twice as fast as the "primesW" version due to the increased algorithmic efficiency!

 

In order to implement a Priority Queue version with Haskell, an efficient Priority Queue, which is not part of the standard Haskell libraries, is required. A Min Heap implementation is likely best suited for this task in providing the most efficient frequently used peeks of the next item in the queue and replacement of the first item in the queue (not using a "pop" followed by a "push) with "pop" operations then not used at all, and "push" operations used relatively infrequently. Judging by O'Neill's use of an efficient "deleteMinAndInsert" operation which she states "(We provide deleteMinAndInsert becausea heap-based implementation can support this operation with considerably less rearrangement than a deleteMin followed by an insert.)", which statement is true for a Min Heap Priority Queue and not others, and her reference to a priority queue by (Paulson, 1996), the queue she used is likely the one as provided as a simple true functional [http://rosettacode.org/wiki/Priority_queue#Haskell Min Heap implementation on RosettaCode], from which the essential functions are reproduced here:

| Br !k v !(PriorityQ k v) !(PriorityQ k v)

deriving (Eq, Ord, Read, Show)

replaceMinPQ :: Ord k => k -> v -> PriorityQ k v -> PriorityQ k v

replaceMinPQ wk wv Mt = Mt

 

The "peekMin" function retrieves both the key and value in a tuple so processing is required to access whichever is required for further processing. As well, the output of the peekMin function is a Maybe with the case of an empty queue providing a Nothing output.

 

The following code is O'Neill's original odds-only code (without wheel factorization) from her paper slightly adjusted as per the requirements of this Min Heap implementation as laid out above; note the `seq` adjustments to the "adjust" function to make the evaluation of the entry tuple more strict for better efficiency:

-- this copyright message is retained and changed versions of the file

-- are clearly marked.

| otherwise = table

where (n, n':ns) = case peekMinPQ table of

 

The above code is almost four times slower than the version of the Tree Merging sieve above for the first million primes although it is about the same speed as the original Richard Bird sieve with the "odds-only" adaptation as above. It is slow and uses a huge amount of memory for primarily one reason: over eagerness in adding prime composite streams to the queue, which are added as the primes are listed rather than when they are required as the output primes stream reaches the square of a given base prime; this over eagerness also means that the processed numbers must have a large range in order to not overflow when squared (as in the default Integer = infinite precision integers as used here and by O'Neill, but Int64's or Word64's would give a practical range) which processing of wide range numbers adds processing and memory requirement overhead. Although O'Neill's code is elegant, it also loses some efficiency due to the extensive use of lazy list processing, not all of which is required even for a wheel factorization implementation.

 

The following code is adjusted to reduce the amount of lazy list processing and to add a secondary base primes stream (or a succession of streams when the combinator is used) so as to overcome the above problems and reduce memory consumption to only that required for the primes below the square root of the currently sieved number; using this means that 32-bit Int's are sufficient for a reasonable range and memory requirements become relatively negligible:

primesPQx() = 2 : _Y ((3 :) . sieve 5 emptyPQ 9) -- initBasePrms

where

in ntbl `seq` adjust ntbl

| otherwise = table

 

The above code is over five times faster than the previous (O'Neill) Priority Queue code half again faster than the Tree-Merging Odds-Only code for a range of a hundred million primes; it is likely faster as the Min Heap is slightly more efficient than Tree Merging due to better tree balancing.

Since the Tree-Folding version above includes the minor changes to work with a factorization wheel, this should have the same minor modifications for comparison purposes, with the code as follows:

 

 

primesPQWheeled :: () -> [Int]

in ntbl `seq` adjust ntbl

| otherwise = table

 

Compiled with -O2 optimization and -fllvm (the LLVM back end), this code gains about the expected ratio in performance in sieving to a range of a hundred million, sieving to this range in about five seconds on a modern medium range desktop computer. This is likely the fastest purely functional incremental type SoE useful for moderate ranges up to about a hundred million to a billion.

 

All of the above unbounded sieves are quite limited in practical sieving range due to the large constant factor overheads in computation, making them mostly just interesting intellectual exercises other than for small ranges of up to about the first million to ten million primes; the following '''"odds-only"''' page-segmented version using (bit-packed internally) mutable unboxed arrays is about 50 times faster than the fastest of the above algorithms for ranges of about that and higher, making it practical for the first several hundred million primes:

 

import Data.Int ( Int64 )

return cmpsts

pagesFrom lwi bps = map (`makePg` bps)

 

The above code as written has a maximum practical range of about 10^12 or so in about an hour.

To show the limitations of the individual prime enumeration, the following code has been refactored from the above to provide an alternate very fast method of counting the unset bits in the culled array (the primes = none composite) using a CPU native pop count instruction:

 

{-# OPTIONS_GHC -O2 -fllvm #-} -- use LLVM for about double speed!

 

putStr $ "Found " ++ show answr

putStr $ " primes up to " ++ show limit

 

When compiled with the "fast way" commented out and the "slow way enabled, the time to find the number of primes up to one billion is about 3.65 seconds on an Intel Sandy Bridge i3-2100 at 3.1 Ghz; with the "fast way" enabled instead, the time is only about 1.45 seconds for the same range, both compiled with the LLVM back end. This shows that more than half of the time for the "slow way" is spent just producing and enumerating the list of primes!

===APL-style===

Rolling set subtraction over the rolling element-wise addition on integers. Basic, slow, worse than quadratic in the number of primes produced, empirically:

. scanl1 minus

Or, a wee bit faster:

. tail) $ a)

A bit optimized, much faster, with better complexity,

. unfoldr (\(a:b:t) -> Just . second ((:t) . (`minus` b))

. span (< head b) $ a)

. scanl1 (zipWith (+) . tail) $ tails [1..]

 

getting nearer to the functional equivalent of the <code>primesPE</code> above, i.e.

. unfoldr (\(a:b:t) -> Just . second ((:t) . (`minus` b))

. span (< head b) $ a)

 

An illustration:

[ 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]

[ 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24, 26, 28, 30, 32]

[ 64, 72, 80, 88, 96,104,112,120,128,136,144,152,160,168,176]

[ 81, 90, 99,108,117,126,135,144,153,162,171,180,189,198,207]

 

=={{header|HicEst}}==

 

sieve = $ > 1 ! = 0 1 1 1 1 ...

DO i = 1, N

IF( sieve(i) ) WRITE() i

 

=={{header|Icon}} and {{header|Unicon}}==

sieve(100)

end

do p[j] := 0

every write(i:=2 to n & p[i] == 1 & i)

 

Alternatively using sets

sieve(100)

end

delete(primes,1)

every write(!sort(primes))

 

=={{header|J}}==

{{eff note|J|i.&.(p:inv) }}

 

end.

 

 

 

 

end.

 

 

 

 

 

 

 

 

=={{header|Java}}==

{{works with|Java|1.5+}}

 

public class Sieve{

return primes;

}

 

To optimize by testing only odd numbers, replace the loop marked "unoptimized" with these lines:

for(int i = 3;i <= n;i += 2){

nums.add(i);

 

Version using a BitSet:

import java.util.BitSet;

 

return primes;

}

 

 

Version using a TreeSet:

import java.util.TreeSet;

 

return list;

}

 

===Infinite iterator===

{{trans|Python}}

{{works with|Java|1.5+}}

import java.util.PriorityQueue;

import java.math.BigInteger;

}

}

{{out}}

<pre>

{{trans|Python}}

{{works with|Java|1.5+}}

import java.util.HashMap;

}

 

{{out}}<pre>Found 5761455 primes up to 100000000 in 4297 milliseconds.</pre>

{{trans|JavaScript}}

{{works with|Java|1.5+}}

import java.util.ArrayList;

 

}

 

{{out}}<pre>Found 50847534 primes up to 1000000000 in 3201 milliseconds.</pre>

 

=={{header|JavaScript}}==

var primes = [];

if (limit >= 2) {

if (typeof print == "undefined")

print = (typeof WScript != "undefined") ? WScript.Echo : alert;

outputs:

<pre>2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97</pre>

Substituting the following code for the function for '''an odds-only algorithm using bit packing''' for the array produces code that is many times faster than the above:

 

var prms = [];

if (limit >= 2) prms = [2];

}

return prms;

 

While the above code is quite fast especially using an efficient JavaScript engine such as Google Chrome's V8, it isn't as elegant as it could be using the features of the new EcmaScript6 specification when it comes out about the end of 2014 and when JavaScript engines including those of browsers implement that standard in that we might choose to implement an incremental algorithm iterators or generators similar to as implemented in Python or F# (yield). Meanwhile, we can emulate some of those features by using a simulation of an iterator class (which is easier than using a call-back function) for an '''"infinite" generator based on an Object dictionary''' as in the following odds-only code written as a JavaScript class:

 

function SoEIncClass() {

this.n = 0;

};

return SoEIncClass;

 

The above code can be used to find the nth prime (which would require estimating the required range limit using the previous fixed range code) by using the following code:

 

for (var i = 1; i < 1000000; i++, gen.next());

var prime = gen.next();

if (typeof print == "undefined")

print = (typeof WScript != "undefined") ? WScript.Echo : alert;

 

to produce the following output (about five seconds using Google Chrome's V8 JavaScript engine):

This can be implemented as '''an "infinite" odds-only generator using page segmentation''' for a considerable speed-up with the alternate JavaScript class code as follows:

 

function SoEPgClass() {

this.bi = -1; // constructor resets the enumeration to start...

};

return SoEPgClass;

 

The above code is about fifty times faster (about five seconds to calculate 50 million primes to about a billion on the Google Chrome V8 JavaScript engine) than the above dictionary based code.

 

The speed for both of these "infinite" solutions will also respond to further wheel factorization techniques, especially for the dictionary based version where any added overhead to deal with the factorization wheel will be negligible compared to the dictionary overhead. The dictionary version would likely speed up about a factor of three or a little more with maximum wheel factorization applied; the page segmented version probably won't gain more than a factor of two and perhaps less due to the overheads of array look-up operations.

 

=={{header|JOVIAL}}==

START

FILE MYOUTPUT ... $ ''Insufficient information to complete this declaration''

END

TERM$

 

=={{header|jq}}==

{{works with|jq|1.4}}

Short and sweet ...

 

# eratosthenes/0 produces an array of primes less than or equal to $n:

def eratosthenes:

 

def erase(i):

else .

end;

| [null, null, range(2; $n)]

| reduce (2, 1 + (2 * range(1; $s))) as $i (.; erase($i))

'''Examples''':

{{out}}

 

[2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97]

{{out}}

664579

 

=={{header|Julia}}==

 

Started with 2 already in the array, and then test only for odd numbers and push the prime ones onto the array.

function sieve(lim :: Int)

if lim < 2 return [] end

end

return result

 

Alternate version using <code>findall</code> to get all primes at once in the end

 

for p in 2:n

end

 

At about 35 seconds for a range of a billion on my Intel Atom i5-Z8350 CPU at 1.92 GHz (single threaded) or about 70 CPU clock cycles per culling operation, the above examples are two of the very slowest ways to compute the Sieve of Eratosthenes over any kind of a reasonable range due to a couple of factors:

If one is going to "crib" the MatLab algorithm as above, one may as well do it using odds-only as per the MatLab built-in. The following alternate code improves on the "Alternate" example above by making it sieve odds-only and adjusting the result array contents after to suit:

 

ni = (n - 1) ÷ 2

isprime = trues(ni)

end

end

 

This takes less about 18.5 seconds or 36 CPU cycles per culling operation to find the primes to a billion, but that is still quite slow compared to what can be done below. Note that the result array needs to be created then copied, created by the <code>findall</code> function, then modified in place by the <code>map!</code> function to transform the indices to primes, and finally copied by the <code>pushfirst!</code> function to add the only even prime of two to the beginning, but these operations are quire fast. However, this still consumes a lot of memory, as in about 64 Megabytes for the sieve buffer and over 400 Megabytes for the result (8-byte Int's for 64 bit execution) to sieve to a billion, and culling the huge culling buffer that doesn't fit the CPU cache sizes is what makes it slow.

The creation of an output results array is not necessary if the purpose is just to scan across the resulting primes once, they can be output using an iterator (from a `BitArray`) as in the following odds-only code:

 

 

struct Primes

return (i + i + 1, i + 1)

end

 

for which using the following code:

 

@time length(Primes(100)) # warm up JIT

# println(@time count(x->true, Primes(1000000000))) # about 1.5 seconds slower counting over iteration

println(@time length(Primes(1000000000)))

end

 

results in the following output:

For any kind of reasonably large range such as a billion, a page segmented version should be used with the pages sized to the CPU caches for much better memory access times. As well, the following odds-only example uses a custom bit packing algorithm for a further two times speed-up, also reducing the memory allocation delays by reusing the sieve buffers when possible (usually possible):

 

const BasePrime = UInt32

const BasePrimesArray = Array{BasePrime,1}

end

low + state + state, state + 1

 

When tested with the following code:

 

print("( ")

for p in PrimesPaged() p > 100 && break; print(p, " ") end

println(@time countPrimesTo(Prime(1000000000)))

end

 

it produces the following:

One of the best simple purely functional Sieve of Eratosthenes algorithms is the infinite tree folding sequence algorithm as implemented in Haskell. As Julia does not have a standard LazyList implementation or library and as a full memoizing lazy list is not required for this algorithm, the following odds-only code implements the rudiments of a Co-Inductive Stream (CIS) in its implementation:

 

 

struct CIS{T}

|> allmultiples |> composites))

CIS{Int}(2, () -> oddprimes())

 

when tested with the following:

 

 

it outputs the following:

To gain some extra speed above the purely functional algorithm above, the Python'ish version as a mutable iterator embedding a mutable standard base Dictionary can be used. The following version uses a secondary delayed injection stream of "base" primes defined recursively to provide the successions of composite values in the Dictionary to be used for sieving:

 

abstract type PrimesDictAbstract end # used for forward reference

mutable struct PrimesDict <: PrimesDictAbstract

end

end

 

The above version can be used and tested with similar code as for the functional version, but is about ten times faster at about 400 CPU clock cycles per culling operation, meaning it has a practical range ten times larger although it still has a O(n (log n) (log log n)) asymptotic performance complexity; for larger ranges such as sieving to a billion or more, this is still over a hundred times slower than the page segmented version using a page segmented sieving array.

 

=={{header|Klingphix}}==

 

%limit %i

pstack

 

 

=={{header|Kotlin}}==

 

fun sieve(max: Int): List<Int> {

fun main(args: Array<String>) {

println(sieve(100))

 

{{out}}

The following code overcomes most of those problems: It only culls odd composites; it culls a bit-packed primitive array (also saving memory); It uses tailcall recursive functions for the loops, which are compiled into simple loops. It also outputs the results as an enumeration, which isn't fast but does not consume any more memory than the culling array. In this way, the program is only limited in sieving range by the maximum size limit of the culling array, although as it grows larger than the CPU cache sizes, it loses greatly in speed; however, that doesn't matter so much if just enumerating the results.

 

val topi = (rng - 3) shr 1

val lstw = topi shr 5

println()

println(primesOdds(1000000).count())

{{output}}

<pre>2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97

Ah, one might say, for such a trivial range one writes for conciseness and not for speed. Well, I say, one can still save memory and some time using odds-only and a bit-packed array, but write very clear and concise (but slower) code using nothing but higher order functions and function calling. The following code using such techniques can use the same "main" function for the same output but is about two times slower, mostly due to the extra time spent making (nested) function calls, including the function calls necessary for enumeration. Note that the effect of using the "(l .. h).forEach { .. }" is the same as the "for i in l .. h { .. }" as both use an iteration across the range but the second is just syntax sugar to make it look more imperative:

 

val topi = (rng - 3) / 2 //convert to nearest index

val size = topi / 32 + 1 //word size to include index

val orng = (-1 .. topi).filter { it < 0 || is_p(it) }.map { i2p(it) }

return Iterable { -> orng.iterator() }

 

The trouble with the above version is that, at least for Kotlin version 1.0, the ".filter" and ".map" extension functions for Iterable<Int> create Java "ArrayList"'s as their output (which are wrapped to return the Kotlin "List<Int>" interface), thus take a considerable amount of memory worse than the first version (using an ArrayList to store the resulting primes), since as the calculations are chained to ".map", require a second ArrayList of up to the same size while the mapping is being done. The following version uses Sequences , which aren't backed by any permanent structure, but it is another small factor slower due to the nested function calls:

 

val topi = (rng - 3) / 2 //convert to nearest index

val size = topi / 32 + 1 //word size to include index

val oseq = iseq(topi, -1).filter { it < 0 || is_p(it) }.map { i2p(it) }

return Iterable { -> oseq.iterator() }

 

===Unbounded Versions===

 

The following Sieve of Eratosthenes is not purely functional in that it uses a Mutable HashMap to store the state of succeeding composite numbers to be skipped over, but embodies the principles of an incremental implementation of the Sieve of Eratosthenes sieving odds-only and is faster than most incremental sieves due to using mutability. As with the fastest of this kind of sieve, it uses a delayed secondary primes feed as a source of base primes to generate the composite number progressions. The code as follows:

yield(2)

fun oddprms(): Sequence<Int> = sequence {

}

yieldAll(oddprms())

 

At about 370 clock cycles per culling operation (about 3,800 cycles per prime) on my tablet class Intel CPU, this is not blazing fast but adequate for ranges of a few millions to a hundred million and thus fine for doing things like solving Euler problems. For instance, Euler Problem 10 of summing the primes to two million can be done with the following "one-liner":

 

to output the correct answer of the following in about 270 milliseconds for my Intel x5-Z8350 at 1.92 Gigahertz:

{{trans|Haskell}}

 

fun toSequence() = generateSequence(this) { it.tailf() } .map { it.head }

}

val stop = System.currentTimeMillis()

println("Took ${stop - strt} milliseconds.")

 

The code is about five times slower than the more imperative hash table based version immediately above due to the costs of the extra levels of function calls in the functional style. The Haskell version from which this is derived is much faster due to the extensive optimizations it does to do with function/closure "lifting" as well as a Garbage Collector specifically tuned for functional code.

 

The very fastest implementations of a primes sieve are all based on bit-packed mutable arrays which can be made unbounded by setting them up so that they are a succession of sieved bit-packed arrays that have been culled of composites. The following code is an odds=only implementation that, again, uses a secondary feed of base primes that is only expanded as necessary (in this case memoized by a rudimentary lazy list structure to avoid recalculation for every base primes sweep per page segment):

internal typealias BasePrime = Int

internal typealias BasePrimeArray = IntArray

}

}

 

For this implementation, counting the primes to a million is trivial at about 15 milliseconds on the same CPU as above, or almost too short to count.

 

Kotlin isn't really fast even as compared to other virtual machine languages such as C# and F# on CLI but that is mostly due to limitations of the Java Virtual Machine (JVM) as to speed of generated Just In Time (JIT) compilation, handling of primitive number operations, enforced array bounds checks, etc. It will always be much slower than native code producing compilers and the (experimental) native compiler for Kotlin still isn't up to speed (pun intended), producing code that is many times slower than the code run on the JVM (December 2018).

 

=={{header|langur}}==

{{trans|D}}

 

 

}

}

}

 

}

 

 

{{out}}

<pre>[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]</pre>

 

=={{header|Liberty BASIC}}==

'The sieve() function uses the flags() array.

'This is a Sieve benchmark adapted from BYTE 1985

end if

next i

 

=={{header|Limbo}}==

 

include "sys.m";

}

}

 

=={{header|Lingo}}==

property _sieve

 

end repeat

end repeat

 

 

{{out}}

 

=={{header|LiveCode}}==

set itemdel to comma

local sieve

sort items of sieve ascending numeric

return sieve

 

 

 

# calculates prime numbers up to a given number

end mouseUp

 

 

[http://www.melellington.com/sieve/livecode-sieve-output.png LiveCode output example]

due to the use of mod (modulo = division) in the filter function.

A coinduction based solution just for fun:

:- object(sieve).

 

 

:- end_object.

Example query:

?- sieve::primes(20, P).

P = [2, 3|_S1], % where

_S1 = [5, 7, 11, 13, 17, 19, 2, 3|_S1] .

 

=={{header|LOLCODE}}==

CAN HAS STDIO?

 

I IZ Eratosumthin YR 100 MKAY

 

 

{{Out}}

 

=={{header|Lua}}==

if n < 2 then return {} end

local t = {0} -- clears '1'

for i = 2, n do if t[i] ~= 0 then table.insert(primes, i) end end

return primes

 

The following changes the code to '''odds-only''' using the same large array-based algorithm:

if n < 2 then return {} end

if n < 3 then return {2} end

for i = 0, lmt do if t[i] ~= 0 then table.insert(primes, i + i + 3) end end

return primes

 

The following code implements '''an odds-only "infinite" generator style using a table as a hash table''', including postponing adding base primes to the table:

 

local _cand = 0; local _lstbp = 3; local _lstsqr = 9

local _composites = {}; local _bps = nil

while gen.next() <= 10000000 do count = count + 1 end -- sieves to 10 million

print(count)

 

which outputs "664579" in about three seconds. As this code uses much less memory for a given range than the previous ones and retains efficiency better with range, it is likely more appropriate for larger sieve ranges.

 

=={{header|M2000 Interpreter}}==

Module EratosthenesSieve (x) {

\\ Κόσκινο του Ερατοσθένη