Sieve of Eratosthenes: Difference between revisions
→{{header|langur}}
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
Langurmonkey (talk | contribs) |
||
| (66 intermediate revisions by 29 users not shown) | |||
Line 761:
ENDLOOP.
</syntaxhighlight>
=={{header|ABC}}==
<syntaxhighlight lang="ABC">HOW TO SIEVE UP TO n:
SHARE sieve
PUT {} IN sieve
FOR cand IN {2..n}: PUT 1 IN sieve[cand]
FOR cand IN {2..floor root n}:
IF sieve[cand] = 1:
PUT cand*cand IN comp
WHILE comp <= n:
PUT 0 IN sieve[comp]
PUT comp+cand IN comp
HOW TO REPORT prime n:
SHARE sieve
IF n<2: FAIL
REPORT sieve[n] = 1
SIEVE UP TO 100
FOR n IN {1..100}:
IF prime n: WRITE n</syntaxhighlight>
{{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|ACL2}}==
Line 1,287 ⟶ 1,310:
The required list of prime divisors obtains by recursion (<tt>{⍵/⍳⍴⍵}∇⌈⍵*0.5</tt>).
=== Optimized Version ===|
<syntaxhighlight lang="apl">sieve←{
Line 2,248 ⟶ 2,271:
@ ^ p3\" ":<
2 234567890123456789012345678901234567890123456789012345678901234567890123456789
=={{header|Binary Lambda Calculus}}==
The BLC sieve of Eratosthenes as documented at https://github.com/tromp/AIT/blob/master/characteristic_sequences/primes.lam is the 167 bit program
<pre>00010001100110010100011010000000010110000010010001010111110111101001000110100001110011010000000000101101110011100111111101111000000001111100110111000000101100000110110</pre>
The infinitely long output is
<pre>001101010001010001010001000001010000010001010001000001000001010000010001010000010001000001000000010001010001010001000000000000010001000001010000000001010000010000010001000001000001010000000001010001010000000000010000000000010001010001000001010000000001000001000001000001010000010001010000000001000000000000010001010001000000000000010000010000000001010001000001000...</pre>
=={{header|BQN}}==
Line 2,765 ⟶ 2,798:
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):
<syntaxhighlight lang="csharp"> static void Main(string[] args) {
Console.WriteLine(new PrimesXXX().ElementAt(1000000 - 1)); // zero based indexing...
}</syntaxhighlight>
Line 2,778 ⟶ 2,811:
<syntaxhighlight lang="cpp">#include <iostream>
#include <iterator>
#include <algorithm>
#include <vector>
// Fills the range [start, end) with 1 if the integer corresponding to the index is composite and 0 otherwise.
// requires: I is RandomAccessIterator
template<typename I>
void mark_composites(I start, I end)
{
std::fill(start, end, 0);
for (
{
if (*
{
auto prime = std::
// mark all multiples of this prime number as composite.
auto multiple_it = it;
while (std::distance(multiple_it, end) > prime)
{
std::advance(multiple_it, prime);
*multiple_it = 1;
}
}
}
}
// Fills "out" with the prime numbers in the range 2...N where N = distance(start, end).
// requires: I is a RandomAccessIterator
// O is an OutputIterator
template <typename I, typename O>
O sieve_primes(I start, I end, O out)
{
mark_composites(start, end);
for (auto it = start + 1; it != end; ++it)
{
if (*it == 0)
{
*
++
}
}
return
}
int main(
{
std::vector<
sieve_primes(is_composite.begin(), is_composite.end(), std::ostream_iterator<int>(std::cout, " "));
// sieve_primes(is_composite.begin(), is_composite.end(), std::back_inserter(primes));
}
=== Boost ===
Line 4,658 ⟶ 4,699:
4987
4999</pre>
=={{header|Craft Basic}}==
<syntaxhighlight lang="basic">define limit = 120
dim flags[limit]
for n = 2 to limit
let flags[n] = 1
next n
print "prime numbers less than or equal to ", limit ," are:"
for n = 2 to sqrt(limit)
if flags[n] = 1 then
for i = n * n to limit step n
let flags[i] = 0
next i
endif
next n
for n = 1 to limit
if flags[n] then
print n
endif
next n</syntaxhighlight>
{{out| Output}}<pre>
prime numbers less than or equal to 120 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 101 103 107 109 113 </pre>
=={{header|Crystal}}==
Line 5,476 ⟶ 5,557:
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|dc}}==
<syntaxhighlight lang="dc">[dn[,]n dsx [d 1 r :a lx + d ln!<.] ds.x lx] ds@
[sn 2 [d;a 0=@ 1 + d ln!<#] ds#x] se
100 lex</syntaxhighlight>
{{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|Delphi}}==
Line 5,706 ⟶ 5,797:
=={{header|EasyLang}}==
<syntaxhighlight lang="text">
len is_divisible[] 100
max = sqrt len is_divisible[]
if
for i =
.
for i = 2 to len is_divisible[]
if is_divisible[i] = 0
.</syntaxhighlight>
Line 6,365 ⟶ 6,451:
=={{header|Elm}}==
===Elm with immutable arrays===
<syntaxhighlight lang="elm">
module PrimeArray exposing (main)
import Array exposing (Array, foldr, map, set)
import Html exposing (div, h1, p, text)
import Html.Attributes exposing (style)
{-
The Eratosthenes sieve task in Rosetta Code does not accept the use of modulo function (allthough Elm functions modBy and remainderBy work always correctly as they require type Int excluding type Float). Thus the solution needs an indexed work array as Elm has no indexes for lists.
In this method we need no division remainder calculations, as we just set the markings of non-primes into the array. We need the indexes that we know, where the marking of the non-primes shall be set.
Because everything is immutable in Elm, every change of array values will create a new array save the original array unchanged. That makes the program running slower or consuming more space of memory than with non-functional imperative languages. All conventional loops (for, while, until) are excluded in Elm because immutability requirement.
Live: https://ellie-app.com/pTHJyqXcHtpa1
-}
alist =
List.range 2 150
-- Work array contains integers 2 ... 149
workArray =
Array.fromList alist
n : Int
n =
List.length alist
-- The max index of integers used in search for primes
-- limit * limit < n
-- equal: limit <= √n
limit : Int
limit =
round (0.5 + sqrt (toFloat n))
-- Remove zero cells of the array
findZero : Int -> Bool
findZero =
\el -> el > 0
zeroFree : Array Int
zeroFree =
Array.filter findZero workResult
nrFoundPrimes =
Array.length zeroFree
workResult : Array Int
workResult =
loopI 2 limit workArray
{- As Elm has no loops (for, while, until)
we must use recursion instead!
The search of prime starts allways saving the
first found value (not setting zero) and continues setting the multiples of prime to zero.
Zero is no integer and may thus be used as marking of non-prime numbers. At the end, only the primes remain in the array and the zeroes are removed from the resulted array to be shown in Html.
-}
-- The recursion increasing variable i follows:
loopI : Int -> Int -> Array Int -> Array Int
loopI i imax arr =
if i > imax then
arr
else
let
arr2 =
phase i arr
in
loopI (i + 1) imax arr2
-- The helper function
phase : Int -> Array Int -> Array Int
phase i =
arrayMarker i (2 * i - 2) n
lastPrime =
Maybe.withDefault 0 <| Array.get (nrFoundPrimes - 1) zeroFree
outputArrayInt : Array Int -> String
outputArrayInt arr =
decorateString <|
foldr (++) "" <|
Array.map (\x -> String.fromInt x ++ " ") arr
decorateString : String -> String
decorateString str =
"[ " ++ str ++ "]"
-- Recursively marking the multiples of p with zero
-- This loop operates with constant p
arrayMarker : Int -> Int -> Int -> Array Int -> Array Int
arrayMarker p min max arr =
let
arr2 =
set min 0 arr
min2 =
min + p
in
if min < max then
arrayMarker p min2 max arr2
else
arr
main =
div [ style "margin" "2%" ]
[ h1 [] [ text "Sieve of Eratosthenes" ]
, text ("List of integers [2, ... ," ++ String.fromInt n ++ "]")
, p [] [ text ("Total integers " ++ String.fromInt n) ]
, p [] [ text ("Max prime of search " ++ String.fromInt limit) ]
, p [] [ text ("The largest found prime " ++ String.fromInt lastPrime) ]
, p [ style "color" "blue", style "font-size" "1.5em" ]
[ text (outputArrayInt zeroFree) ]
, p [] [ text ("Found " ++ String.fromInt nrFoundPrimes ++ " primes") ]
]
</syntaxhighlight>
{{out}}
<pre>
List of integers [2, ... ,149]
Total integers 149
Max prime of search 13
The largest found prime 149
[ 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 101 103 107 109 113 127 131 137 139 149 ]
Found 35 primes </pre>
===Concise Elm Immutable Array Version===
Although functional, the above code is written in quite an imperative style, so the following code is written in a more concise functional style and includes timing information for counting the number of primes to a million:
<syntaxhighlight lang="elm">module Main exposing (main)
import Browser exposing (element)
import Task exposing (Task, succeed, perform, andThen)
import Html exposing (Html, div, text)
import Time exposing (now, posixToMillis)
import Array exposing (repeat, get, set)
cLIMIT : Int
cLIMIT = 1000000
primesArray : Int -> List Int
primesArray n =
if n < 2 then [] else
let
sz = n + 1
loopbp bp arr =
let s = bp * bp in
if s >= sz then arr else
let tst = get bp arr |> Maybe.withDefault True in
if tst then loopbp (bp + 1) arr else
let
cullc c iarr =
if c >= sz then iarr else
cullc (c + bp) (set c True iarr)
in loopbp (bp + 1) (cullc s arr)
cmpsts = loopbp 2 (repeat sz False)
cnvt (i, t) = if t then Nothing else Just i
in cmpsts |> Array.toIndexedList
|> List.drop 2 -- skip the values for zero and one
|> List.filterMap cnvt -- primes are indexes of not composites
type alias Model = List String
type alias Msg = Model
test : (Int -> List Int) -> Int -> Cmd Msg
test primesf lmt =
let
to100 = primesf 100 |> List.map String.fromInt |> String.join ", "
to100str = "The primes to 100 are: " ++ to100
timemillis() = now |> andThen (succeed << posixToMillis)
in timemillis() |> andThen (\ strt ->
let cnt = primesf lmt |> List.length
in timemillis() |> andThen (\ stop ->
let answrstr = "Found " ++ (String.fromInt cnt) ++ " primes to "
++ (String.fromInt cLIMIT) ++ " in "
++ (String.fromInt (stop - strt)) ++ " milliseconds."
in succeed [to100str, answrstr] ) ) |> perform identity
main : Program () Model Msg
main =
element { init = \ _ -> ( [], test primesArray cLIMIT )
, update = \ msg _ -> (msg, Cmd.none)
, subscriptions = \ _ -> Sub.none
, view = div [] << List.map (div [] << List.singleton << text) }</syntaxhighlight>
{{out}}
<pre>The primes up to 100 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.
Found 78498 primes to 1000000 in 958 milliseconds.</pre>
The above output is the contents of the HTML web page as shown with Google Chrome version 1.23 running on an AMD 7840HS CPU at 5.1 GHz (single thread boosted).
===Concise Elm Immutable Array Odds-Only Version===
The following code can replace the `primesArray` function in the above program and called from the testing and display code (two places):
<syntaxhighlight lang="elm">primesArrayOdds : Int -> List Int
primesArrayOdds n =
if n < 2 then [] else
let
sz = (n - 1) // 2
loopi i arr =
let s = (i + i) * (i + 3) + 3 in
if s >= sz then arr else
let tst = get i arr |> Maybe.withDefault True in
if tst then loopi (i + 1) arr else
let
bp = i + i + 3
cullc c iarr =
if c >= sz then iarr else
cullc (c + bp) (set c True iarr)
in loopi (i + 1) (cullc s arr)
cmpsts = loopi 0 (repeat sz False)
cnvt (i, t) = if t then Nothing else Just <| i + i + 3
oddprms = cmpsts |> Array.toIndexedList |> List.filterMap cnvt
in 2 :: oddprms</syntaxhighlight>
{{out}}
<pre>The primes up to 100 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.
Found 78498 primes to 1000000 in 371 milliseconds.</pre>
The above output is the contents of the HTML web page as shown with Google Chrome version 1.23 running on an AMD 7840HS CPU at 5.1 GHz (single thread boosted).
===Richard Bird Tree Folding Elm Version===
The Elm language doesn't efficiently handle the Sieve of Eratosthenes (SoE) algorithm because it doesn't have directly accessible linear arrays (the Array module used above is based on a persistent tree of sub arrays) and also does Copy On Write (COW) for every write to every location as well as a logarithmic process of updating as a "tree" to minimize the COW operations. Thus, there is better performance implementing the Richard Bird Tree Folding functional algorithm, as follows:
{{trans|Haskell}}
<syntaxhighlight lang="elm">module Main exposing (main)
Line 6,373 ⟶ 6,725:
import Html exposing (Html, div, text)
import Time exposing (now, posixToMillis)
cLIMIT : Int
cLIMIT = 1000000
type CIS a = CIS a (() -> CIS a)
Line 6,415 ⟶ 6,770:
in CIS 2 <| \ () -> oddprms()
type alias Model =
type alias Msg =
test : (() -> CIS Int) -> Int -> Cmd Msg
test primesf lmt =
let
to100 = primesf() |> uptoCIS2List 100
|> List.map String.fromInt |> String.join ", "
to100str = "The primes to 100 are: " ++ to100
timemillis() = now |> andThen (succeed << posixToMillis)
in timemillis() |> andThen (\ strt ->
let cnt =
in timemillis() |> andThen (\
let answrstr = "Found " ++ (String.fromInt cnt) ++ " primes to "
++ (String.fromInt cLIMIT) ++ " in "
++ (String.fromInt (stop - strt)) ++ " milliseconds."
in succeed [to100str, answrstr] ) ) |> perform identity
main : Program () Model Msg
main =
element { init = \ _ -> (
, update = \ msg _ ->
, subscriptions = \ _ -> Sub.none
, view =
{{out}}
<pre>The primes up to 100 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.
Found 78498 primes to 1000000 in
The above output is the contents of the HTML web page as shown with Google Chrome version 1.23 running on an AMD 7840HS CPU at 5.1 GHz (single thread boosted).
===Elm Priority Queue Version===
Using a Binary Minimum Heap Priority Queue is a constant factor faster than the above code as the data structure is balanced rather than "heavy to the right" and requires less memory allocations/deallocation in the following code, which implements enough of the Priority Queue for the purpose. Just substitute the following code for `primesTreeFolding` and pass `primesPQ` as an argument to `test` rather than `primesTreeFolding`:
<syntaxhighlight lang="elm">type PriorityQ comparable v =
Mt
| Br comparable v (PriorityQ comparable v)
Line 6,549 ⟶ 6,863:
else CIS n <| \ () -> sieve (n + 2) pq q bps
oddprms() = CIS 3 <| \ () -> sieve 5 emptyPQ 9 <| oddprms()
in CIS 2 <| \ () -> oddprms()</syntaxhighlight>
{{out}}
<pre>The primes up to 100 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.
Found 78498 primes to 1000000 in
The above output is the contents of the HTML web page as shown with Google Chrome version 1.23 running on an AMD 7840HS CPU at 5.1 GHz (single thread boosted).
=={{header|Emacs Lisp}}==
Line 6,887 ⟶ 7,163:
16301 16319 16333 16339 16349 16361 16363 16369 16381
1899 PRIMES
</pre>
=={{header|Euler}}==
The original Euler doesn't have loops built-in. Loops can easily be added by defining and calling suitable procedures with literal procedures as parameters. In this sample, a C-style "for" loop procedure is defined and used to sieve and print the primes.<br>
'''begin'''
'''new''' sieve; '''new''' for; '''new''' prime; '''new''' i;
for <- ` '''formal''' init; '''formal''' test; '''formal''' incr; '''formal''' body;
'''begin'''
'''label''' again;
init;
again: '''if''' test '''then''' '''begin''' body; incr; '''goto''' again '''end''' '''else''' 0
'''end'''
'
;
sieve <- ` '''formal''' n;
'''begin'''
'''new''' primes; '''new''' i; '''new''' i2; '''new''' j;
primes <- '''list''' n;
for( ` i <- 1 ', ` i <= n ', ` i <- i + 1 '
, ` primes[ i ] <- '''true''' '
);
primes[ 1 ] <- '''false''';
for( ` i <- 2 '
, ` [ i2 <- i * i ] <= n '
, ` i <- i + 1 '
, ` '''if''' primes[ i ] '''then'''
for( ` j <- i2 ', ` j <= n ', ` j <- j + i '
, ` primes[ j ] <- '''false''' '
)
'''else''' 0
'
);
primes
'''end'''
'
;
prime <- sieve( 30 );
for( ` i <- 1 ', ` i <= '''length''' prime ', ` i <- i + 1 '
, ` '''if''' prime[ i ] '''then''' '''out''' i '''else''' 0 '
)
'''end''' $
{{out}}
<pre>
NUMBER 2
NUMBER 3
NUMBER 5
NUMBER 7
NUMBER 11
NUMBER 13
NUMBER 17
NUMBER 19
NUMBER 23
NUMBER 29
</pre>
Line 8,156 ⟶ 8,490:
I: DWord;
begin
Write(aList[I], ', ');
WriteLn(aList[aList.Size - 1]);
end;
aList.Free;
end;
Line 8,297 ⟶ 8,634:
for I := 0 to Length(aList) - 2 do
Write(aList[I], ', ');
WriteLn(aList[High(aList)]);
end;
begin
Line 8,406 ⟶ 8,744:
</syntaxhighlight>
=={{header|Peri}}==
''Note: With benchmark function''
<syntaxhighlight lang="peri">
###sysinclude standard.uh
tick sto startingtick
#g 100000 sto MAX
@MAX mem !maximize sto primeNumbers
one count
2 0 sto#s primeNumbers
2 @MAX külső: {{ ,
@count {{
{{}}§külső primeNumbers[{{}}] !/ else {{<}}§külső
}} // @count vége
//{{}} gprintnl // A talált prímszám kiiratásához kommentezzük ki e sort
{{}} @count++ sto#s primeNumbers
}} // @MAX vége
@primeNumbers inv mem
//."Time : " tick @startingtick - print ." tick\n"
."Prímek száma = " @count printnl
end
{ „MAX” } { „startingtick” } { „primeNumbers” } { „count” }
</syntaxhighlight>
=={{header|FutureBasic}}==
Line 8,438 ⟶ 8,799:
=={{header|Fōrmulæ}}==
{{FormulaeEntry|page=https://formulae.org/?script=examples/Sieve_of_Eratosthenes}}
'''Solution'''
[[File:Fōrmulæ - Sieve of Eratosthenes 01.png]]
'''Test case'''
[[File:Fōrmulæ - Sieve of Eratosthenes 02.png]]
[[File:Fōrmulæ - Sieve of Eratosthenes 03.png]]
=={{header|GAP}}==
Line 8,841 ⟶ 9,208:
var p int
p = <-primes
Line 9,062 ⟶ 9,429:
a = accumArray (\ b c -> False) True (1,q-1)
[(i,()) | (s,y) <- fs, i <- [y+s, y+s+s..q]]</syntaxhighlight>
====As list comprehension====
<syntaxhighlight lang="haskell">import Data.Array.Unboxed
import Data.List (tails, inits)
primes = 2 : [ n |
(r:q:_, px) <- zip (tails (2 : [p*p | p <- primes]))
(inits primes),
(n, True) <- assocs ( accumArray (\_ _ -> False) True
(r+1,q-1)
[ (m,()) | p <- px
, s <- [ div (r+p) p * p]
, m <- [s,s+p..q-1] ] :: UArray Int Bool
) ]</syntaxhighlight>
===Basic list-based sieve===
Line 10,301 ⟶ 10,683:
=={{header|jq}}==
{{works with|jq|1.4}}
==Bare Bones==
Short and sweet ...
Line 10,308 ⟶ 10,690:
def eratosthenes:
def erase(i):
if .[i] then
reduce (range(2*i; length; i)) as $j (.; .[$j] = false)
else .
end;
Line 10,327 ⟶ 10,709:
{{out}}
664579
===Enhanced Sieve===
Here is a more economical variant that:
* produces a stream of primes less than or equal to a given integer;
* only records the status of odd integers greater than 3 during the sieving process;
* optimizes the inner loop as described in the task description.
<syntaxhighlight lang=jq>
def primes:
# The array we use for the sieve only stores information for the odd integers greater than 1:
# index integer
# 0 3
# k 2*k + 3
# So if we wish to mark m = 2*k + 3, the relevant index is: m - 3 / 2
def ix:
if . % 2 == 0 then null
else ((. - 3) / 2)
end;
# erase(i) sets .[i*j] to false for odd integral j > i, and assumes i is odd
def erase(i):
((i - 3) / 2) as $k
# Consider relevant multiples:
then (((length * 2 + 3) / i)) as $upper
# ... only consider odd multiples from i onwards
| reduce range(i; $upper; 2) as $j (.;
(((i * $j) - 3) / 2) as $m
| if .[$m] then .[$m] = false else . end);
if . < 2 then []
else (. + 1) as $n
| (($n|sqrt) / 2) as $s
| [range(3; $n; 2)|true]
| reduce (1 + (2 * range(1; $s)) ) as $i (.; erase($i))
| . as $sieve
| 2, (range(3; $n; 2) | select($sieve[ix]))
end ;
def count(s): reduce s as $_ (0; .+1);
count(1e6 | primes)
</syntaxhighlight>
{{output}}
<pre>
78498
</pre>
=={{header|Julia}}==
Line 10,354 ⟶ 10,783:
Alternate version using <code>findall</code> to get all primes at once in the end
<syntaxhighlight lang="julia">function sieve(n
for p in 2:n
primes[p .* (2:n÷p)]
end
findall(primes)
end</syntaxhighlight>
Line 11,228 ⟶ 11,650:
=={{header|Lambdatalk}}==
<syntaxhighlight lang="scheme">
• 1) create an array of natural numbers, [0,1,2,3, ... ,n-1]
• 2) the 3rd number is 2, we set to dots all its composites by steps of 2,
• 3) the 4th number is 3, we set to dots all its composites by steps of 3,
• 4) the 6th number is 5, we set to dots all its composites by steps of 5,
• 5) the remaining numbers are primes and we clean all dots.
For instance:
1: 0 0 0 0 0 0 0 0 0 9 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 3
0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0
2: 0 1 2 3 . 5 . 7 . 9 . 1 . 3 . 5 . 7 . 9 . 1 . 3 . 5 . 7 . 9 .
3: 0 1 2 3 . 5 . 7 . . . 1 . 3 . . . 7 . 9 . . . 3 . 5 . . . 9 .
4: 0 1 2 3 . 5 . 7 . .
| | | | | | | | | |
5: 0 0 0 0 1 1 1 1 2 2
2 3 5 7 1 3 7 9 3 9
1) recursive version {rsieve n}
{def rsieve
{def rsieve.mark
{lambda {:n :a :i :j}
{if {< :j :n}
then {rsieve.mark :n
{A.set! :j . :a}
:i
{+ :i :j}}
else :a}}}
{def rsieve.loop
{lambda {:n :a :i}
{if {< {* :i :i} :n}
then {rsieve.loop :n
{if {W.equal? {A.get :i :a} .}
then :a
else {rsieve.mark :n :a :i {* :i :i}}}
{+ :i 1}}
else :a}}}
{lambda {:n}
{S.replace \s by space in
{S.replace (\[|\]|\.|,) by space in
{A.disp
{A.slice 2 -1
{rsieve.loop :n
{A.new {S.serie 0 :n}} 2}}}}}}}
-> rsieve
{rsieve 1000}
-> 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 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 223 227 229 233 239 241 251 257 263 269 271 277 281 283 293 307 311 313 317 331 337 347 349 353 359 367 373 379 383 389 397 401 409 419 421 431 433 439 443 449 457 461 463 467 479 487 491 499 503 509 521 523 541 547 557 563 569 571 577 587 593 599 601 607 613 617 619 631 641 643 647 653 659 661 673 677 683 691 701 709 719 727 733 739 743 751 757 761 769 773 787 797 809 811 821 823 827 829 839 853 857 859 863 877 881 883 887 907 911 919 929 937 941 947 953 967 971 977 983 991 997
Note: this version doesn't avoid stackoverflow.
2) iterative version {isieve n}
{def isieve
{def isieve.mark
{lambda {:n :a :i}
{S.map {{lambda {:a :j} {A.set! :j . :a}
} :a}
{S.serie {* :i :i} :n :i} }}}
{lambda {:n}
{S.replace \s by space in
{S.replace (\[|\]|\.|,) by space in
{A.disp
{A.slice 2 -1
{S.last
{S.map {{lambda {:n :a :i} {if {W.equal? {A.get :i :a} .}
then
else {isieve.mark :n :a :i}}
} :n {A.new {S.serie 0 :n}}}
{S.serie 2 {sqrt :n} 1}}}}}}}}}
-> isieve
{isieve 1000}
-> 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 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 223 227 229 233 239 241 251 257 263 269 271 277 281 283 293 307 311 313 317 331 337 347 349 353 359 367 373 379 383 389 397 401 409 419 421 431 433 439 443 449 457 461 463 467 479 487 491 499 503 509 521 523 541 547 557 563 569 571 577 587 593 599 601 607 613 617 619 631 641 643 647 653 659 661 673 677 683 691 701 709 719 727 733 739 743 751 757 761 769 773 787 797 809 811 821 823 827 829 839 853 857 859 863 877 881 883 887 907 911 919 929 937 941 947 953 967 971 977 983 991 997
Notes:
- this version avoids stackoverflow.
- From 1 to 1000000 there are 78500 primes (computed in ~15000ms) and the last is 999983.
</syntaxhighlight>
=={{header|langur}}==
{{trans|D}}
<syntaxhighlight lang="langur">val .sieve =
if .limit < 2:
var .composite =
.composite[1] = true
for .n in 2
if not .composite[.n] {
for .k = .n^2 ; .k < .limit ; .k += .n {
.composite[.k] = true
}
Line 11,279 ⟶ 11,754:
}
}
writeln .sieve(100)
</syntaxhighlight>
{{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|LFE}}==
<syntaxhighlight lang="lisp">
(defmodule eratosthenes
(export (sieve 1)))
(defun sieve (limit)
(sieve limit (lists:seq 2 limit)))
(defun sieve
((limit (= l (cons p _))) (when (> (* p p) limit))
l)
((limit (cons p ns))
(cons p (sieve limit (remove-multiples p (* p p) ns)))))
(defun remove-multiples (p q l)
(lists:reverse (remove-multiples p q l '())))
(defun remove-multiples
((_ _ '() s) s)
((p q (cons q ns) s)
(remove-multiples p q ns s))
((p q (= r (cons a _)) s) (when (> a q))
(remove-multiples p (+ q p) r s))
((p q (cons n ns) s)
(remove-multiples p q ns (cons n s))))
</syntaxhighlight>
{{Out}}
<pre>
lfe> (slurp "sieve.lfe")
#(ok eratosthenes)
lfe> (sieve 100)
(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}}==
Line 11,688 ⟶ 12,198:
Module EratosthenesSieve (x) {
\\ Κόσκινο του Ερατοσθένη
Dim i(x+1): k=2: k2=sqrt(x)
While k<=k2{if i(k) else for m=k*k to x step k{i(m)=1}
Print str$(timecount/1000,"####0.##")+" s"
Input "Press enter skip print or a non zero to get results:", a%
if a% then For i=2to
}
Print:Print "Done"
}
EratosthenesSieve 1000
Line 12,130 ⟶ 12,634:
IO.Put("\n");
END Sieve.</syntaxhighlight>
=={{header|Mojo}}==
Tested with Mojo version 0.7:
<syntaxhighlight lang="mojo">from memory import memset_zero
from memory.unsafe import (DTypePointer)
from time import (now)
alias cLIMIT: Int = 1_000_000_000
struct SoEBasic(Sized):
var len: Int
var cmpsts: DTypePointer[DType.bool] # because DynamicVector has deep copy bug in mojo version 0.7
var sz: Int
var ndx: Int
fn __init__(inout self, limit: Int):
self.len = limit - 1
self.sz = limit - 1
self.ndx = 0
self.cmpsts = DTypePointer[DType.bool].alloc(limit - 1)
memset_zero(self.cmpsts, limit - 1)
for i in range(limit - 1):
let s = i * (i + 4) + 2
if s >= limit - 1: break
if self.cmpsts[i]: continue
let bp = i + 2
for c in range(s, limit - 1, bp):
self.cmpsts[c] = True
for i in range(limit - 1):
if self.cmpsts[i]: self.sz -= 1
fn __del__(owned self):
self.cmpsts.free()
fn __copyinit__(inout self, existing: Self):
self.len = existing.len
self.cmpsts = DTypePointer[DType.bool].alloc(self.len)
for i in range(self.len):
self.cmpsts[i] = existing.cmpsts[i]
self.sz = existing.sz
self.ndx = existing.ndx
fn __moveinit__(inout self, owned existing: Self):
self.len = existing.len
self.cmpsts = existing.cmpsts
self.sz = existing.sz
self.ndx = existing.ndx
fn __len__(self: Self) -> Int: return self.sz
fn __iter__(self: Self) -> Self: return self
fn __next__(inout self: Self) -> Int:
if self.ndx >= self.len: return 0
while (self.ndx < self.len) and (self.cmpsts[self.ndx]):
self.ndx += 1
let rslt = self.ndx + 2; self.sz -= 1; self.ndx += 1
return rslt
fn main():
print("The primes to 100 are:")
for prm in SoEBasic(100): print_no_newline(prm, " ")
print()
let strt0 = now()
let answr0 = len(SoEBasic(1_000_000))
let elpsd0 = (now() - strt0) / 1000000
print("Found", answr0, "primes up to 1,000,000 in", elpsd0, "milliseconds.")
let strt1 = now()
let answr1 = len(SoEBasic(cLIMIT))
let elpsd1 = (now() - strt1) / 1000000
print("Found", answr1, "primes up to", cLIMIT, "in", elpsd1, "milliseconds.")</syntaxhighlight>
{{out}}
<pre>The primes to 100 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
Found 78498 primes up to 1,000,000 in 1.2642770000000001 milliseconds.
Found 50847534 primes up to 1000000000 in 6034.328751 milliseconds.</pre>
as run on an AMD 7840HS CPU at 5.1 GHz.
Note that due to the huge memory array used, when large ranges are selected, the speed is disproportional in speed slow down by about four times.
This solution uses an interator struct which seems to be the Mojo-preferred way to do this, and normally a DynamicVector would have been used the the culling array except that there is a bug in this version of DynamicVector where the array is not properly deep copied when copied to a new location, so the raw pointer type is used.
===Odds-Only with Optimizations===
This version does three significant improvements to the above code as follows:
1) It is trivial to skip the processing to store representations for and cull the even comosite numbers other than the prime number two, saving half the storage space and reducing the culling time to about 40 percent.
2) There is a repeating pattern of culling composite representations over a bit-packed byte array (which reduces the storage requirement by another eight times) that repeats every eight culling operations, which can be encapsulated by a extreme loop unrolling technique with compiler generated constants as done here.
3) Further, there is a further extreme optimization technique of dense culling for small base prime values whose culling span is less than one register in size where the loaded register is repeatedly culled for different base prime strides before being written out (with such optimization done by the compiler), again using compiler generated modification constants. This technique is usually further optimizated by modern compilers to use efficient autovectorization and the use of SIMD registers available to the architecture to reduce these culling operations to an avererage of a tiny fraction of a CPU clock cycle per cull.
Mojo version 0.7 was tested:
<syntaxhighlight lang="mojo">from memory import (memset_zero, memcpy)
from memory.unsafe import (DTypePointer)
from math.bit import ctpop
from time import (now)
alias cLIMIT: Int = 1_000_000_000
alias cBufferSize: Int = 262144 # bytes
alias cBufferBits: Int = cBufferSize * 8
alias UnrollFunc = fn(DTypePointer[DType.uint8], Int, Int, Int) -> None
@always_inline
fn extreme[OFST: Int, BP: Int](pcmps: DTypePointer[DType.uint8], bufsz: Int, s: Int, bp: Int):
var cp = pcmps + (s >> 3)
let r1: Int = ((s + bp) >> 3) - (s >> 3)
let r2: Int = ((s + 2 * bp) >> 3) - (s >> 3)
let r3: Int = ((s + 3 * bp) >> 3) - (s >> 3)
let r4: Int = ((s + 4 * bp) >> 3) - (s >> 3)
let r5: Int = ((s + 5 * bp) >> 3) - (s >> 3)
let r6: Int = ((s + 6 * bp) >> 3) - (s >> 3)
let r7: Int = ((s + 7 * bp) >> 3) - (s >> 3)
let plmt: DTypePointer[DType.uint8] = pcmps + bufsz - r7
while cp < plmt:
cp.store(cp.load() | (1 << OFST))
(cp + r1).store((cp + r1).load() | (1 << ((OFST + BP) & 7)))
(cp + r2).store((cp + r2).load() | (1 << ((OFST + 2 * BP) & 7)))
(cp + r3).store((cp + r3).load() | (1 << ((OFST + 3 * BP) & 7)))
(cp + r4).store((cp + r4).load() | (1 << ((OFST + 4 * BP) & 7)))
(cp + r5).store((cp + r5).load() | (1 << ((OFST + 5 * BP) & 7)))
(cp + r6).store((cp + r6).load() | (1 << ((OFST + 6 * BP) & 7)))
(cp + r7).store((cp + r7).load() | (1 << ((OFST + 7 * BP) & 7)))
cp += bp
let eplmt: DTypePointer[DType.uint8] = plmt + r7
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << OFST))
cp += r1
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + BP) & 7)))
cp += r2 - r1
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 2 * BP) & 7)))
cp += r3 - r2
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 3 * BP) & 7)))
cp += r4 - r3
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 4 * BP) & 7)))
cp += r5 - r4
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 5 * BP) & 7)))
cp += r6 - r5
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 6 * BP) & 7)))
cp += r7 - r6
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 7 * BP) & 7)))
fn mkExtrm[CNT: Int](pntr: Pointer[UnrollFunc]):
@parameter
if CNT >= 32:
return
alias OFST = CNT >> 2
alias BP = ((CNT & 3) << 1) + 1
pntr.offset(CNT).store(extreme[OFST, BP])
mkExtrm[CNT + 1](pntr)
@always_inline
fn mkExtremeFuncs() -> Pointer[UnrollFunc]:
let jmptbl: Pointer[UnrollFunc] = Pointer[UnrollFunc].alloc(32)
mkExtrm[0](jmptbl)
return jmptbl
let extremeFuncs = mkExtremeFuncs()
alias DenseFunc = fn(DTypePointer[DType.uint64], Int, Int) -> DTypePointer[DType.uint64]
fn mkDenseCull[N: Int, BP: Int](cp: DTypePointer[DType.uint64]):
@parameter
if N >= 64:
return
alias MUL = N * BP
var cop = cp.offset(MUL >> 6)
cop.store(cop.load() | (1 << (MUL & 63)))
mkDenseCull[N + 1, BP](cp)
@always_inline
fn denseCullFunc[BP: Int](pcmps: DTypePointer[DType.uint64], bufsz: Int, s: Int) -> DTypePointer[DType.uint64]:
var cp: DTypePointer[DType.uint64] = pcmps + (s >> 6)
let plmt = pcmps + (bufsz >> 3) - BP
while cp < plmt:
mkDenseCull[0, BP](cp)
cp += BP
return cp
fn mkDenseFunc[CNT: Int](pntr: Pointer[DenseFunc]):
@parameter
if CNT >= 64:
return
alias BP = (CNT << 1) + 3
pntr.offset(CNT).store(denseCullFunc[BP])
mkDenseFunc[CNT + 1](pntr)
@always_inline
fn mkDenseFuncs() -> Pointer[DenseFunc]:
let jmptbl : Pointer[DenseFunc] = Pointer[DenseFunc].alloc(64)
mkDenseFunc[0](jmptbl)
return jmptbl
let denseFuncs : Pointer[DenseFunc] = mkDenseFuncs()
@always_inline
fn cullPass(cmpsts: DTypePointer[DType.uint8], bytesz: Int, s: Int, bp: Int):
if bp <= 129: # dense culling
var sm = s
while (sm >> 3) < bytesz and (sm & 63) != 0:
cmpsts[sm >> 3] |= (1 << (sm & 7))
sm += bp
let bcp = denseFuncs[(bp - 3) >> 1](cmpsts.bitcast[DType.uint64](), bytesz, sm)
var ns = 0
var ncp = bcp
let cmpstslmtp = (cmpsts + bytesz).bitcast[DType.uint64]()
while ncp < cmpstslmtp:
ncp[0] |= (1 << (ns & 63))
ns += bp
ncp = bcp + (ns >> 6)
else: # extreme loop unrolling culling
extremeFuncs[((s & 7) << 2) + ((bp & 7) >> 1)](cmpsts, bytesz, s, bp)
# for c in range(s, self.len, bp): # slow bit twiddling way
# self.cmpsts[c >> 3] |= (1 << (c & 7))
fn countPagePrimes(ptr: DTypePointer[DType.uint8], bitsz: Int) -> Int:
let wordsz: Int = (bitsz + 63) // 64 # round up to nearest 64 bit boundary
var rslt: Int = wordsz * 64
let bigcmps = ptr.bitcast[DType.uint64]()
for i in range(wordsz - 1):
rslt -= ctpop(bigcmps[i]).to_int()
rslt -= ctpop(bigcmps[wordsz - 1] | (-2 << ((bitsz - 1) & 63))).to_int()
return rslt
struct SoEOdds(Sized):
var len: Int
var cmpsts: DTypePointer[DType.uint8] # because DynamicVector has deep copy bug in Mojo version 0.7
var sz: Int
var ndx: Int
fn __init__(inout self, limit: Int):
self.len = 0 if limit < 2 else (limit - 3) // 2 + 1
self.sz = 0 if limit < 2 else self.len + 1 # for the unprocessed only even prime of two
self.ndx = -1
let bytesz = 0 if limit < 2 else ((self.len + 63) & -64) >> 3 # round up to nearest 64 bit boundary
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
memset_zero(self.cmpsts, bytesz)
for i in range(self.len):
let s = (i + i) * (i + 3) + 3
if s >= self.len: break
if (self.cmpsts[i >> 3] >> (i & 7)) & 1 != 0: continue
let bp = i + i + 3
cullPass(self.cmpsts, bytesz, s, bp)
self.sz = countPagePrimes(self.cmpsts, self.len) + 1 # add one for only even prime of two
fn __del__(owned self):
self.cmpsts.free()
fn __copyinit__(inout self, existing: Self):
self.len = existing.len
let bytesz = (self.len + 7) // 8
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
memcpy(self.cmpsts, existing.cmpsts, bytesz)
self.sz = existing.sz
self.ndx = existing.ndx
fn __moveinit__(inout self, owned existing: Self):
self.len = existing.len
self.cmpsts = existing.cmpsts
self.sz = existing.sz
self.ndx = existing.ndx
fn __len__(self: Self) -> Int: return self.sz
fn __iter__(self: Self) -> Self: return self
@always_inline
fn __next__(inout self: Self) -> Int:
if self.ndx < 0:
self.ndx = 0; self.sz -= 1; return 2
while (self.ndx < self.len) and ((self.cmpsts[self.ndx >> 3] >> (self.ndx & 7)) & 1 != 0):
self.ndx += 1
let rslt = (self.ndx << 1) + 3; self.sz -= 1; self.ndx += 1; return rslt
fn main():
print("The primes to 100 are:")
for prm in SoEOdds(100): print_no_newline(prm, " ")
print()
let strt0 = now()
let answr0 = len(SoEOdds(1_000_000))
let elpsd0 = (now() - strt0) / 1000000
print("Found", answr0, "primes up to 1,000,000 in", elpsd0, "milliseconds.")
let strt1 = now()
let answr1 = len(SoEOdds(cLIMIT))
let elpsd1 = (now() - strt1) / 1000000
print("Found", answr1, "primes up to", cLIMIT, "in", elpsd1, "milliseconds.")</syntaxhighlight>
{{out}}
<pre>The primes to 100 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
Found 78498 primes up to 1,000,000 in 0.085067000000000004 milliseconds.
Found 50847534 primes up to 1000000000 in 1204.866606 milliseconds.</pre>
This was run on the same computer as the above example; notice that while this is much faster than that version, it is still very slow as the sieving range gets large such that the relative processing time for a range that is 1000 times as large is about ten times slower than as might be expected by simple scaling. This is due to the "one huge sieving buffer" algorithm that gets very large with increasing range (and in fact will eventually limit the sieving range that can be used) to exceed the size of CPU cache buffers and thus greatly slow average memory access times.
===Page-Segmented Odds-Only with Optimizations===
While the above version performs reasonably well for small sieving ranges that fit within the CPU caches of a few tens of millions, as one can see it gets much slower with larger ranges and as well its huge RAM memory consumption limits the maximum range over which it can be used. This version solves these problems be breaking the huge sieving array into "pages" that each fit within the CPU cache size and processing each "page" sequentially until the target range is reached. This technique also greatly reduces memory requirements to only that required to store the base prime value representations up to the square root of the range limit (about O(n/log n) storage plus a fixed size page buffer. In this case, the storage for the base primes has been reduced by a constant factor by storing them as single byte deltas from the previous value, which works for ranges up to the 64-bit number range where the biggest gap is two times 192 and since we store only for odd base primes, the gap values are all half values to fit in a single byte.
Currently, Mojo has problems with some functions in the standard libraries such as the integer square root function is not accurate nor does it work for the required integer types so a custom integer square root function is supplied. As well, current Mojo does not support recursion for hardly any useful cases (other than compile time global function recursion), so the `SoeOdds` structure from the previous answer had to be kept to generate the base prime representation table (or this would have had to be generated from scratch within the new `SoEOddsPaged` structure). Finally, it didn't seem to be worth using the `Sized` trait for the new structure as this would seem to sometimes require processing the pages twice, one to obtain the size and once if iteration across the prime values is required.
Tested with Mojo version 0.7:
<syntaxhighlight lang="mojo">from memory import (memset_zero, memcpy)
from memory.unsafe import (DTypePointer)
from math.bit import ctpop
from time import (now)
alias cLIMIT: Int = 1_000_000_000
alias cBufferSize: Int = 262144 # bytes
alias cBufferBits: Int = cBufferSize * 8
fn intsqrt(n: UInt64) -> UInt64:
if n < 4:
if n < 1: return 0 else: return 1
var x: UInt64 = n; var qn: UInt64 = 0; var r: UInt64 = 0
while qn < 64 and (1 << qn) <= n:
qn += 2
var q: UInt64 = 1 << qn
while q > 1:
if qn >= 64:
q = 1 << (qn - 2); qn = 0
else:
q >>= 2
let t: UInt64 = r + q
r >>= 1
if x >= t:
x -= t; r += q
return r
alias UnrollFunc = fn(DTypePointer[DType.uint8], Int, Int, Int) -> None
@always_inline
fn extreme[OFST: Int, BP: Int](pcmps: DTypePointer[DType.uint8], bufsz: Int, s: Int, bp: Int):
var cp = pcmps + (s >> 3)
let r1: Int = ((s + bp) >> 3) - (s >> 3)
let r2: Int = ((s + 2 * bp) >> 3) - (s >> 3)
let r3: Int = ((s + 3 * bp) >> 3) - (s >> 3)
let r4: Int = ((s + 4 * bp) >> 3) - (s >> 3)
let r5: Int = ((s + 5 * bp) >> 3) - (s >> 3)
let r6: Int = ((s + 6 * bp) >> 3) - (s >> 3)
let r7: Int = ((s + 7 * bp) >> 3) - (s >> 3)
let plmt: DTypePointer[DType.uint8] = pcmps + bufsz - r7
while cp < plmt:
cp.store(cp.load() | (1 << OFST))
(cp + r1).store((cp + r1).load() | (1 << ((OFST + BP) & 7)))
(cp + r2).store((cp + r2).load() | (1 << ((OFST + 2 * BP) & 7)))
(cp + r3).store((cp + r3).load() | (1 << ((OFST + 3 * BP) & 7)))
(cp + r4).store((cp + r4).load() | (1 << ((OFST + 4 * BP) & 7)))
(cp + r5).store((cp + r5).load() | (1 << ((OFST + 5 * BP) & 7)))
(cp + r6).store((cp + r6).load() | (1 << ((OFST + 6 * BP) & 7)))
(cp + r7).store((cp + r7).load() | (1 << ((OFST + 7 * BP) & 7)))
cp += bp
let eplmt: DTypePointer[DType.uint8] = plmt + r7
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << OFST))
cp += r1
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + BP) & 7)))
cp += r2 - r1
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 2 * BP) & 7)))
cp += r3 - r2
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 3 * BP) & 7)))
cp += r4 - r3
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 4 * BP) & 7)))
cp += r5 - r4
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 5 * BP) & 7)))
cp += r6 - r5
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 6 * BP) & 7)))
cp += r7 - r6
if eplmt == cp or eplmt < cp: return
cp.store(cp.load() | (1 << ((OFST + 7 * BP) & 7)))
fn mkExtrm[CNT: Int](pntr: Pointer[UnrollFunc]):
@parameter
if CNT >= 32:
return
alias OFST = CNT >> 2
alias BP = ((CNT & 3) << 1) + 1
pntr.offset(CNT).store(extreme[OFST, BP])
mkExtrm[CNT + 1](pntr)
@always_inline
fn mkExtremeFuncs() -> Pointer[UnrollFunc]:
let jmptbl: Pointer[UnrollFunc] = Pointer[UnrollFunc].alloc(32)
mkExtrm[0](jmptbl)
return jmptbl
let extremeFuncs = mkExtremeFuncs()
alias DenseFunc = fn(DTypePointer[DType.uint64], Int, Int) -> DTypePointer[DType.uint64]
fn mkDenseCull[N: Int, BP: Int](cp: DTypePointer[DType.uint64]):
@parameter
if N >= 64:
return
alias MUL = N * BP
var cop = cp.offset(MUL >> 6)
cop.store(cop.load() | (1 << (MUL & 63)))
mkDenseCull[N + 1, BP](cp)
@always_inline
fn denseCullFunc[BP: Int](pcmps: DTypePointer[DType.uint64], bufsz: Int, s: Int) -> DTypePointer[DType.uint64]:
var cp: DTypePointer[DType.uint64] = pcmps + (s >> 6)
let plmt = pcmps + (bufsz >> 3) - BP
while cp < plmt:
mkDenseCull[0, BP](cp)
cp += BP
return cp
fn mkDenseFunc[CNT: Int](pntr: Pointer[DenseFunc]):
@parameter
if CNT >= 64:
return
alias BP = (CNT << 1) + 3
pntr.offset(CNT).store(denseCullFunc[BP])
mkDenseFunc[CNT + 1](pntr)
@always_inline
fn mkDenseFuncs() -> Pointer[DenseFunc]:
let jmptbl : Pointer[DenseFunc] = Pointer[DenseFunc].alloc(64)
mkDenseFunc[0](jmptbl)
return jmptbl
let denseFuncs : Pointer[DenseFunc] = mkDenseFuncs()
@always_inline
fn cullPass(cmpsts: DTypePointer[DType.uint8], bytesz: Int, s: Int, bp: Int):
if bp <= 129: # dense culling
var sm = s
while (sm >> 3) < bytesz and (sm & 63) != 0:
cmpsts[sm >> 3] |= (1 << (sm & 7))
sm += bp
let bcp = denseFuncs[(bp - 3) >> 1](cmpsts.bitcast[DType.uint64](), bytesz, sm)
var ns = 0
var ncp = bcp
let cmpstslmtp = (cmpsts + bytesz).bitcast[DType.uint64]()
while ncp < cmpstslmtp:
ncp[0] |= (1 << (ns & 63))
ns += bp
ncp = bcp + (ns >> 6)
else: # extreme loop unrolling culling
extremeFuncs[((s & 7) << 2) + ((bp & 7) >> 1)](cmpsts, bytesz, s, bp)
# for c in range(s, self.len, bp): # slow bit twiddling way
# self.cmpsts[c >> 3] |= (1 << (c & 7))
fn cullPage(lwi: Int, lmt: Int, cmpsts: DTypePointer[DType.uint8], bsprmrps: DTypePointer[DType.uint8]):
var bp = 1; var ndx = 0
while True:
bp += bsprmrps[ndx].to_int() << 1
let i = (bp - 3) >> 1
var s = (i + i) * (i + 3) + 3
if s >= lmt: break
if s >= lwi: s -= lwi
else:
s = (lwi - s) % bp
if s != 0: s = bp - s
cullPass(cmpsts, cBufferSize, s, bp)
ndx += 1
fn countPagePrimes(ptr: DTypePointer[DType.uint8], bitsz: Int) -> Int:
let wordsz: Int = (bitsz + 63) // 64 # round up to nearest 64 bit boundary
var rslt: Int = wordsz * 64
let bigcmps = ptr.bitcast[DType.uint64]()
for i in range(wordsz - 1):
rslt -= ctpop(bigcmps[i]).to_int()
rslt -= ctpop(bigcmps[wordsz - 1] | (-2 << ((bitsz - 1) & 63))).to_int()
return rslt
struct SoEOdds(Sized):
var len: Int
var cmpsts: DTypePointer[DType.uint8] # because DynamicVector has deep copy bug in Mojo version 0.7
var sz: Int
var ndx: Int
fn __init__(inout self, limit: Int):
self.len = 0 if limit < 2 else (limit - 3) // 2 + 1
self.sz = 0 if limit < 2 else self.len + 1 # for the unprocessed only even prime of two
self.ndx = -1
let bytesz = 0 if limit < 2 else ((self.len + 63) & -64) >> 3 # round up to nearest 64 bit boundary
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
memset_zero(self.cmpsts, bytesz)
for i in range(self.len):
let s = (i + i) * (i + 3) + 3
if s >= self.len: break
if (self.cmpsts[i >> 3] >> (i & 7)) & 1 != 0: continue
let bp = i + i + 3
cullPass(self.cmpsts, bytesz, s, bp)
self.sz = countPagePrimes(self.cmpsts, self.len) + 1 # add one for only even prime of two
fn __del__(owned self):
self.cmpsts.free()
fn __copyinit__(inout self, existing: Self):
self.len = existing.len
let bytesz = (self.len + 7) // 8
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
memcpy(self.cmpsts, existing.cmpsts, bytesz)
self.sz = existing.sz
self.ndx = existing.ndx
fn __moveinit__(inout self, owned existing: Self):
self.len = existing.len
self.cmpsts = existing.cmpsts
self.sz = existing.sz
self.ndx = existing.ndx
fn __len__(self: Self) -> Int: return self.sz
fn __iter__(self: Self) -> Self: return self
@always_inline
fn __next__(inout self: Self) -> Int:
if self.ndx < 0:
self.ndx = 0; self.sz -= 1; return 2
while (self.ndx < self.len) and ((self.cmpsts[self.ndx >> 3] >> (self.ndx & 7)) & 1 != 0):
self.ndx += 1
let rslt = (self.ndx << 1) + 3; self.sz -= 1; self.ndx += 1; return rslt
struct SoEOddsPaged:
var len: Int
var cmpsts: DTypePointer[DType.uint8] # because DynamicVector has deep copy bug in Mojo version 0.7
var sz: Int # 0 means finished; otherwise contains number of odd base primes
var ndx: Int
var lwi: Int
var bsprmrps: DTypePointer[DType.uint8] # contains deltas between odd base primes starting from zero
fn __init__(inout self, limit: UInt64):
self.len = 0 if limit < 2 else ((limit - 3) // 2 + 1).to_int()
self.sz = 0 if limit < 2 else 1 # means iterate until this is set to zero
self.ndx = -1 # for unprocessed only even prime of two
self.lwi = 0
if self.len < cBufferBits:
let bytesz = ((self.len + 63) & -64) >> 3 # round up to nearest 64 bit boundary
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
self.bsprmrps = DTypePointer[DType.uint8].alloc(self.sz)
else:
self.cmpsts = DTypePointer[DType.uint8].alloc(cBufferSize)
let bsprmitr = SoEOdds(intsqrt(limit).to_int())
self.sz = len(bsprmitr)
self.bsprmrps = DTypePointer[DType.uint8].alloc(self.sz)
var ndx = -1; var oldbp = 1
for bsprm in bsprmitr:
if ndx < 0: ndx += 1; continue # skip over the 2 prime
self.bsprmrps[ndx] = (bsprm - oldbp) >> 1
oldbp = bsprm; ndx += 1
self.bsprmrps[ndx] = 255 # one extra value to go beyond the necessary cull space
fn __del__(owned self):
self.cmpsts.free(); self.bsprmrps.free()
fn __copyinit__(inout self, existing: Self):
self.len = existing.len
self.sz = existing.sz
let bytesz = cBufferSize if self.len >= cBufferBits
else ((self.len + 63) & -64) >> 3 # round up to nearest 64 bit boundary
self.cmpsts = DTypePointer[DType.uint8].alloc(bytesz)
memcpy(self.cmpsts, existing.cmpsts, bytesz)
self.ndx = existing.ndx
self.lwi = existing.lwi
self.bsprmrps = DTypePointer[DType.uint8].alloc(self.sz)
memcpy(self.bsprmrps, existing.bsprmrps, self.sz)
fn __moveinit__(inout self, owned existing: Self):
self.len = existing.len
self.cmpsts = existing.cmpsts
self.sz = existing.sz
self.ndx = existing.ndx
self.lwi = existing.lwi
self.bsprmrps = existing.bsprmrps
fn countPrimes(self) -> Int:
if self.len <= cBufferBits: return len(SoEOdds(2 * self.len + 1))
var cnt = 1; var lwi = 0
let cmpsts = DTypePointer[DType.uint8].alloc(cBufferSize)
memset_zero(cmpsts, cBufferSize)
cullPage(0, cBufferBits, cmpsts, self.bsprmrps)
while lwi + cBufferBits <= self.len:
cnt += countPagePrimes(cmpsts, cBufferBits)
lwi += cBufferBits
memset_zero(cmpsts, cBufferSize)
let lmt = lwi + cBufferBits if lwi + cBufferBits <= self.len else self.len
cullPage(lwi, lmt, cmpsts, self.bsprmrps)
cnt += countPagePrimes(cmpsts, self.len - lwi)
return cnt
fn __len__(self: Self) -> Int: return self.sz
fn __iter__(self: Self) -> Self: return self
@always_inline
fn __next__(inout self: Self) -> Int: # don't count number of primes by interating - slooow
if self.ndx < 0:
self.ndx = 0; self.lwi = 0
if self.len < 2: self.sz = 0
elif self.len <= cBufferBits:
let bytesz = ((self.len + 63) & -64) >> 3 # round up to nearest 64 bit boundary
memset_zero(self.cmpsts, bytesz)
for i in range(self.len):
let s = (i + i) * (i + 3) + 3
if s >= self.len: break
if (self.cmpsts[i >> 3] >> (i & 7)) & 1 != 0: continue
let bp = i + i + 3
cullPass(self.cmpsts, bytesz, s, bp)
else:
memset_zero(self.cmpsts, cBufferSize)
cullPage(0, cBufferBits, self.cmpsts, self.bsprmrps)
return 2
let rslt = ((self.lwi + self.ndx) << 1) + 3; self.ndx += 1
if self.lwi + cBufferBits >= self.len:
while (self.lwi + self.ndx < self.len) and ((self.cmpsts[self.ndx >> 3] >> (self.ndx & 7)) & 1 != 0):
self.ndx += 1
else:
while (self.ndx < cBufferBits) and ((self.cmpsts[self.ndx >> 3] >> (self.ndx & 7)) & 1 != 0):
self.ndx += 1
while (self.ndx >= cBufferBits) and (self.lwi + cBufferBits <= self.len):
self.ndx = 0; self.lwi += cBufferBits; memset_zero(self.cmpsts, cBufferSize)
let lmt = self.lwi + cBufferBits if self.lwi + cBufferBits <= self.len else self.len
cullPage(self.lwi, lmt, self.cmpsts, self.bsprmrps)
let buflmt = cBufferBits if self.lwi + cBufferBits <= self.len else self.len - self.lwi
while (self.ndx < buflmt) and ((self.cmpsts[self.ndx >> 3] >> (self.ndx & 7)) & 1 != 0):
self.ndx += 1
if self.lwi + self.ndx >= self.len: self.sz = 0
return rslt
fn main():
print("The primes to 100 are:")
for prm in SoEOddsPaged(100): print_no_newline(prm, " ")
print()
let strt0 = now()
let answr0 = SoEOddsPaged(1_000_000).countPrimes()
let elpsd0 = (now() - strt0) / 1000000
print("Found", answr0, "primes up to 1,000,000 in", elpsd0, "milliseconds.")
let strt1 = now()
let answr1 = SoEOddsPaged(cLIMIT).countPrimes()
let elpsd1 = (now() - strt1) / 1000000
print("Found", answr1, "primes up to", cLIMIT, "in", elpsd1, "milliseconds.")</syntaxhighlight>
{{out}}
<pre>The primes to 100 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
Found 78498 primes up to 1,000,000 in 0.084122000000000002 milliseconds.
Found 50847534 primes up to 1000000000 in 139.509275 milliseconds.</pre>
This was tested on the same computer as the previous Mojo versions. Note that the time now scales quite well with range since there are no longer the huge RAM access time bottleneck's. This version is only about 2.25 times slower than Kim Walich's primesieve program written in C++ and the mostly constant factor difference will be made up if one adds wheel factorization to the same level as he uses (basic wheel factorization ratio of 48/105 plus some other more minor optimizations). This version can count the number of primes to 1e11 in about 21.85 seconds on this machine. It will work reasonably efficiently up to a range of about 1e14 before other optimization techniques such as "bucket sieving" should be used.
For counting the number of primes to a billion (1e9), this version has reduced the time by about a factor of 40 from the original version and over eight times from the odds-only version above. Adding wheel factorization will make it almost two and a half times faster yet for a gain in speed of about a hundred times over the original version.
=={{header|MUMPS}}==
Line 13,145 ⟶ 14,286:
<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 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 223 227 229 233 239 241 251 257 263 269 271 277 281 283 293 307 311 313 317 331 337 347 349 353 359 367 373 379 383 389 397 401 409 419 421 431 433 439 443 449 457 461 463 467 479 487 491 499 503 509 521 523 541 547 557 563 569 571 577 587 593 599 601 607 613 617 619 631 641 643 647 653 659 661 673 677 683 691 701 709 719 727 733 739 743 751 757 761 769 773 787 797 809 811 821 823 827 829 839 853 857 859 863 877 881 883 887 907 911 919 929 937 941 947 953 967 971 977 983 991 997)
</pre>
=={{header|ooRexx}}==
<syntaxhighlight lang="oorexx">
/*ooRexx program generates & displays primes via the sieve of Eratosthenes.
* derived from first Rexx version
* uses an array rather than a stem for the list
* uses string methods rather than BIFs
* uses new ooRexx keyword LOOP, extended assignment
* and line comments
* uses meaningful variable names and restructures code
* layout for improved understandability
****************************************************************************/
arg highest --get highest number to use.
if \highest~datatype('W') then
highest = 200 --use default value.
isPrime = .array~new(highest) --container for all numbers.
isPrime~fill(1) --assume all numbers are prime.
w = highest~length --width of the biggest number,
-- it's used for aligned output.
out1 = 'prime'~right(20) --first part of output messages.
np = 0 --no primes so far.
loop j = 2 for highest - 1 --all numbers up through highest.
if isPrime[j] = 1 then do --found one.
np += 1 --bump the prime counter.
say out1 np~right(w) ' --> ' j~right(w) --display output.
loop m = j * j to highest by j
isPrime[m] = '' --strike all multiples: not prime.
end
end
end
say
say np~right(out1~length + 1 + w) 'primes found up to and including ' highest
exit
</syntaxhighlight>
{{out}}
<pre>
prime 1 --> 2
prime 2 --> 3
prime 3 --> 5
prime 4 --> 7
prime 5 --> 11
prime 6 --> 13
prime 7 --> 17
prime 8 --> 19
prime 9 --> 23
prime 10 --> 29
prime 11 --> 31
prime 12 --> 37
prime 13 --> 41
prime 14 --> 43
prime 15 --> 47
prime 16 --> 53
prime 17 --> 59
prime 18 --> 61
prime 19 --> 67
prime 20 --> 71
prime 21 --> 73
prime 22 --> 79
prime 23 --> 83
prime 24 --> 89
prime 25 --> 97
prime 26 --> 101
prime 27 --> 103
prime 28 --> 107
prime 29 --> 109
prime 30 --> 113
prime 31 --> 127
prime 32 --> 131
prime 33 --> 137
prime 34 --> 139
prime 35 --> 149
prime 36 --> 151
prime 37 --> 157
prime 38 --> 163
prime 39 --> 167
prime 40 --> 173
prime 41 --> 179
prime 42 --> 181
prime 43 --> 191
prime 44 --> 193
prime 45 --> 197
prime 46 --> 199
46 primes found up to and including 200
</pre>
===Wheel Version===
<syntaxhighlight lang="oorexx">
/*ooRexx program generates primes via sieve of Eratosthenes algorithm.
* wheel version, 2 handled as special case
* loops optimized: outer loop stops at the square root of
* the limit, inner loop starts at the square of the
* prime just found
* use a list rather than an array and remove composites
* rather than just mark them
* convert list of primes to a list of output messages and
* display them with one say statement
*******************************************************************************/
arg highest -- get highest number to use.
if \highest~datatype('W') then
highest = 200 -- use default value.
w = highest~length -- width of the biggest number,
-- it's used for aligned output.
thePrimes = .list~of(2) -- the first prime is 2.
loop j = 3 to highest by 2 -- populate the list with odd nums.
thePrimes~append(j)
end
j = 3 -- first prime (other than 2)
ix = thePrimes~index(j) -- get the index of 3 in the list.
loop while j*j <= highest -- strike multiples of odd ints.
-- up to sqrt(highest).
loop jm = j*j to highest by j+j -- start at J squared, incr. by 2*J.
thePrimes~removeItem(jm) -- delete it since it's composite.
end
ix = thePrimes~next(ix) -- the index of the next prime.
j = thePrimes[ix] -- the next prime.
end
np = thePrimes~items -- the number of primes since the
-- list is now only primes.
out1 = ' prime number' -- first part of output messages.
out2 = ' --> ' -- middle part of output messages.
ix = thePrimes~first
loop n = 1 to np -- change the list of primes
-- to output messages.
thePrimes[ix] = out1 n~right(w) out2 thePrimes[ix]~right(w)
ix = thePrimes~next(ix)
end
last = np~right(out1~length+1+w) 'primes found up to and including ' highest
thePrimes~append(.endofline || last) -- add blank line and summary line.
say thePrimes~makearray~toString -- display the output.
exit
</syntaxhighlight>
{{out}}when using the limit of 100
<pre>
prime number 1 --> 2
prime number 2 --> 3
prime number 3 --> 5
prime number 4 --> 7
prime number 5 --> 11
prime number 6 --> 13
prime number 7 --> 17
prime number 8 --> 19
prime number 9 --> 23
prime number 10 --> 29
prime number 11 --> 31
prime number 12 --> 37
prime number 13 --> 41
prime number 14 --> 43
prime number 15 --> 47
prime number 16 --> 53
prime number 17 --> 59
prime number 18 --> 61
prime number 19 --> 67
prime number 20 --> 71
prime number 21 --> 73
prime number 22 --> 79
prime number 23 --> 83
prime number 24 --> 89
prime number 25 --> 97
25 primes found up to and including 100
</pre>
Line 15,057 ⟶ 16,360:
number on the primes ancillary stack is set.
Initially all the bits are set except for 0 and 1,
which are not prime numbers by definition.
"eratosthenes" unsets all bits above those specified
by it's argument. )
[ bit ~
Line 15,064 ⟶ 16,369:
[ bit primes share & 0 != ] is isprime ( n --> b )
[ dup dup sqrt times
[ i^ 1+
dup isprime if
[ dup 2 **
[ dup -composite
over +
rot 2dup >
dip unrot until ]
drop ]
drop ]
drop
1+ bit 1 -
primes take &
primes put ] is eratosthenes ( n --> )
100 eratosthenes
Line 15,451 ⟶ 16,759:
== [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]
</syntaxhighlight>
=={{header|Refal}}==
<syntaxhighlight lang="refal">$ENTRY Go {
= <Print <Primes 100>>;
};
Primes {
s.N = <Sieve <Iota 2 s.N>>;
};
Iota {
s.End s.End = s.End;
s.Start s.End = s.Start <Iota <+ 1 s.Start> s.End>;
};
Cross {
s.Step e.List = <Cross (s.Step 1) s.Step e.List>;
(s.Step s.Skip) = ;
(s.Step 1) s.Item e.List = X <Cross (s.Step s.Step) e.List>;
(s.Step s.N) s.Item e.List = s.Item <Cross (s.Step <- s.N 1>) e.List>;
};
Sieve {
= ;
X e.List = <Sieve e.List>;
s.N e.List = s.N <Sieve <Cross s.N e.List>>;
};</syntaxhighlight>
{{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|REXX}}==
Line 15,696 ⟶ 17,033:
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|RPL}}==
This is a direct translation from Wikipedia. The variable <code>i</code> has been renamed <code>ii</code> to avoid confusion with the language constant <code>i</code>=√ -1
{{works with|Halcyon Calc|4.2.8}}
{| class="wikitable"
! RPL code
! Comment
|-
|
≪ → n
≪ { } n + 1 CON 'A' STO
2 n √ '''FOR''' ii
'''IF''' A ii GET '''THEN'''
ii SQ n '''FOR''' j
'A' j 0 PUT ii '''STEP'''
'''END'''
'''NEXT'''
{ }
2 n '''FOR''' ii '''IF''' A ii GET '''THEN''' ii + '''END NEXT'''
'A' PURGE
≫ ≫ '<span style="color:blue">SIEVE</span>' STO
|
<span style="color:blue">SIEVE</span> ''( n -- { prime_numbers } )''
let A be an array of Boolean values, indexed by 2 to n,
initially all set to true.
for i = 2, 3, 4, ..., not exceeding √n do
if A[i] is true
for j = i^2, i^2+i,... not exceeding n do
set A[j] := false
return all i such that A[i] is true.
|}
100 <span style="color:blue">SIEVE</span>
{{out}}
<pre>
1: { 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>
Latest RPL versions allow to remove some slow <code>FOR..NEXT</code> loops and use local variables only.
{{works with|HP|49}}
« 'X' DUP 1 4 PICK 1 SEQ DUP → n a seq123
« 2 n √ '''FOR''' ii
'''IF''' a ii GET '''THEN'''
ii SQ n '''FOR''' j
'a' j 0 PUT ii '''STEP'''
'''END'''
'''NEXT'''
a seq123 IFT TAIL
» » '<span style="color:blue">SIEVE</span>' STO
{{works with|HP|49}}
=={{header|Ruby}}==
Line 16,978 ⟶ 18,369:
Original source: [http://seed7.sourceforge.net/algorith/math.htm#sieve_of_eratosthenes]
=={{header|SETL}}==
<syntaxhighlight lang="setl">program eratosthenes;
print(sieve 100);
op sieve(n);
numbers := [1..n];
numbers(1) := om;
loop for i in [2..floor sqrt n] do
loop for j in [i*i, i*i+i..n] do
numbers(j) := om;
end loop;
end loop;
return [n : n in numbers | n /= om];
end op;
end program;</syntaxhighlight>
{{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|Sidef}}==
Line 17,265 ⟶ 18,674:
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|SparForte}}==
As a structured script.
<syntaxhighlight lang="ada">#!/usr/local/bin/spar
pragma annotate( summary, "sieve" );
pragma annotate( description, "The Sieve of Eratosthenes is a simple algorithm that" );
pragma annotate( description, "finds the prime numbers up to a given integer. Implement ");
pragma annotate( description, "this algorithm, with the only allowed optimization that" );
pragma annotate( description, "the outer loop can stop at the square root of the limit," );
pragma annotate( description, "and the inner loop may start at the square of the prime" );
pragma annotate( description, "just found. That means especially that you shouldn't" );
pragma annotate( description, "optimize by using pre-computed wheels, i.e. don't assume" );
pragma annotate( description, "you need only to cross out odd numbers (wheel based on" );
pragma annotate( description, "2), numbers equal to 1 or 5 modulo 6 (wheel based on 2" );
pragma annotate( description, "and 3), or similar wheels based on low primes." );
pragma annotate( see_also, "http://rosettacode.org/wiki/Sieve_of_Eratosthenes" );
pragma annotate( author, "Ken O. Burtch" );
pragma license( unrestricted );
pragma restriction( no_external_commands );
procedure sieve is
last_bool : constant positive := 20;
type bool_array is array(2..last_bool) of boolean;
a : bool_array;
test_num : positive;
-- limit : positive := positive(numerics.sqrt(float(arrays.last(a))));
-- n : positive := 2;
begin
for i in arrays.first(a)..last_bool loop
a(i) := true;
end loop;
for num in arrays.first(a)..last_bool loop
if a(num) then
test_num := num * num;
while test_num <= last_bool loop
a(test_num) := false;
test_num := @ + num;
end loop;
end if;
end loop;
for i in arrays.first(a)..last_bool loop
if a(i) then
put_line(i);
end if;
end loop;
end sieve;</syntaxhighlight>
=={{header|Standard ML}}==
Line 18,501 ⟶ 19,961:
End Module</syntaxhighlight>{{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 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 211 223 227 229 233 239 241 251 257 263 269 271 277 281 283 293 307 311 313 317 331 337 347 349 353 359 367 373 379 383 389 397 401 409 419 421 431 433 439 443 449 457 461 463 467 479 487 491 499 </pre>
=={{header|V (Vlang)}}==
{{trans|go}}
===Basic sieve of array of booleans===
<syntaxhighlight lang="v (vlang)">fn main() {
limit := 201 // means sieve numbers < 201
Line 18,781 ⟶ 20,241:
=={{header|Wren}}==
<syntaxhighlight lang="
if (n < 2) return []
var comp = List.filled(n-1, false)
| |||