List of primes

primes.md

Primes Benchmark

I got nerd-sniped from this LinkedIn post. On my computer, the R screenshot above ran in 1.564 minutes on my laptop for primes up to 1 million (Macbook Air M1 with arm64 binaries for Julia/R).

I had to give it a shot for comparison's sake. A basic 1:1 translation in Julia was 13.43 seconds (7x faster), and following the same conceptual algorithm (check for remainders up to sqrt(x)) but without as many intermediate allocations was 0.017 seconds (over 5,250x faster).

Prime divisibility check algorithim

R - 93.8 seconds

library(magrittr)
library(purrr)
library(dplyr)
library(lubridate)

primes_less_than <- 1000000

v.primes <- c(2,3)

start <- now()
walk((2:ceiling(primes_less_than/2)) * 2 -1,
  function(x) {
    if ((map_dbl(v.primes[v.primes <= sqrt(x)],
                 ~ x %% .) %>% min) > 0) {
                   v.primes[length(v.primes) + 1] <<- x
    }
    }) %>%
  suppressWarnings

print(now() - start)

Julia - 13.43 sec (7x faster)

function primes_up_to(k) 
	primes = [2,3]
	for n in 3:2:k 
		remainders = n .% primes[primes .<= √(n)]
		if  minimum(remainders;init=1) > 0 
			push!(primes,n)
		end
	end
	return primes
end

Julia - 0.029 seconds (3,250x faster)

This follows the same idea (check for remainders up to sqrt(k)) but does not have as many intermediate allocations.

function primes_up_to(k)
	primes = [2,3]
	for n in 3:2:k 
		d = false
		i = 1
		while (~d) && (i < length(primes)) && (primes[i] <= sqrt(n)) 
			
			if iszero(n % primes[i])
				d = true
			end
			i += 1
		end
		if ~d
			push!(primes,n)
		end
			
			
	end
	return primes
end

Julia - 0.0175 seconds (5,500x faster)

By splitting into sub-functions, it becomes more readable (IMO) and facilitates some further optimizations.

function primes_up_to(k)
	primes = [2,3]
	for n in 3:2:k 
		if ~is_divisible(primes,n)
			push!(primes,n)
		end
	end
	return primes
end

function is_divisible(prime_list,n)
	# integer sqrt keeps everything integers instead of comparing float with int
	limit = isqrt(n) 
	
	for i ∈ eachindex(prime_list)
		p = prime_list[i]
		
		# short circuit and return if the conditionals are met
		p > limit && return false
		iszero(n % p) && return true
					
	end
	
	# we've exhausted the search and can return false
	return false	
end

Other methods

These aren't directly comparable, as they have different algorithms. Two other options:

Primes.jl - 0.0013 seconds

using Primes
primes(1_000_000)

FastPrimeSieve - 0.0003 seconds

using FastPrimeSieve
collect(FastPrimeSieve.SmallSieve(1_000_000))

@bkamins your version is indeed faster:

Main solution:

julia> @benchmark primes_up_to(1_000_000)
BenchmarkTools.Trial: 321 samples with 1 evaluation.
 Range (min … max):  15.367 ms …  21.083 ms  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     15.521 ms               ┊ GC (median):    0.00%
 Time  (mean ± σ):   15.582 ms ± 394.369 μs  ┊ GC (mean ± σ):  0.05% ± 0.56%

    ▁▅▇██▆▆▄▃▂▁
  ▄████████████▆▇▆▄▆▆▄▄▄▁▁▁▁▁▁▄▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▄▁▁▁▁▁▁▄▄▁▁▁▁▄ ▇
  15.4 ms       Histogram: log(frequency) by time      16.8 ms <

 Memory estimate: 1.83 MiB, allocs estimate: 11.

UInt32 version:

julia> @benchmark primes_up_to(1_000_000)
BenchmarkTools.Trial: 361 samples with 1 evaluation.
 Range (min … max):  13.486 ms …  15.491 ms  ┊ GC (min … max): 0.00% … 11.32%
 Time  (median):     13.738 ms               ┊ GC (median):    0.00%
 Time  (mean ± σ):   13.851 ms ± 364.746 μs  ┊ GC (mean ± σ):  0.65% ±  2.20%

      ▁▄▇▇█▆▅▃▁▁
  ▄▄▇▇██████████▇▆▄▁▄▆▆▆▁▄▄▁▄▁▁▁▆▁▁▄▁▁▄▄▄▄▆▆▁▄▁▄▆▁▄▇▁▆▄▄▄▁▄▆▄▆ ▇
  13.5 ms       Histogram: log(frequency) by time      15.2 ms <

 Memory estimate: 3.42 MiB, allocs estimate: 15.