Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!
Get Started for Free

primes

Find all Primes Less Than n


Description

Find all prime numbers aka ‘primes’ less than n.

Uses an obvious sieve method (and some care), working with logical and and integers to be quite fast.

Usage

primes(n, pSeq = NULL)

Arguments

n

a (typically positive integer) number.

pSeq

optionally a vector of primes (2,3,5,...) as if from a primes() call; must be correct. The goal is a speedup, but currently we have not found one single case, where using a non-NULL pSeq is faster.

Details

As the function only uses max(n), n can also be a vector of numbers.

The famous prime number theorem states that π(n), the number of primes below n is asymptotically n / \log(n) in the sense that lim[n -> Inf] π(n) * log(n) / n ~ 1.

Equivalently, the inverse of pi(), the n-th prime number p_n is around n \log n; recent results (Pierre Dusart, 1999), prove that

log n + log log n - 1 < p_n / n < log n + log log n for n >= 6.

Value

numeric vector of all prime numbers <= n.

Author(s)

Bill Venables (<= 2001); Martin Maechler gained another 40% speed, carefully working with logicals and integers.

See Also

factorize. For large n, use the gmp package and its isprime and nextprime functions.

Examples

(p1 <- primes(100))
 system.time(p1k <- primes(1000)) # still lightning fast
 stopifnot(length(p1k) == 168)

 system.time(p.e7 <- primes(1e7)) # still only 0.3 sec (2015 (i7))
 stopifnot(length(p.e7) == 664579)
 ## The famous  pi(n) :=  number of primes <= n:
 pi.n <- approxfun(p.e7, seq_along(p.e7), method = "constant")
 pi.n(c(10, 100, 1000)) # 4 25 168
 plot(pi.n, 2, 1e7, n = 1024, log="xy", axes = FALSE,
      xlab = "n", ylab = quote(pi(n)),
      main = quote("The prime number function " ~ pi(n)))
 eaxis(1); eaxis(2)


## Exploring  p(n) := the n-th prime number  ~=~ n * pnn(n), where
## pnn(n) := log n + log log n
pnn <- function(n) { L <- log(n); L + log(L) }
n <- 6:(N <- length(PR <- primes(1e5)))
m.pn <- cbind(l.pn = ceiling(n*(pnn(n)-1)), pn = PR[n], u.pn = floor(n*pnn(n)))
matplot(n, m.pn, type="l", ylab = quote(p[n]), main = quote(p[n] ~~
        "with lower/upper bounds" ~ n*(log(n) + log(log(n)) -(1~"or"~0))))
## (difference to the lower approximation) / n   --> ~ 0.0426  (?) :
plot(n, PR[n]/n - (pnn(n)-1), type = 'l', cex = 1/8, log="x", xaxt="n")
eaxis(1); abline(h=0, col=adjustcolor(1, 0.5))

sfsmisc

Utilities from 'Seminar fuer Statistik' ETH Zurich

v1.1-11
GPL (>= 2)
Authors
Martin Maechler [aut, cre] (<https://orcid.org/0000-0002-8685-9910>), Werner Stahel [ctb] (Functions: compresid2way(), f.robftest(), last(), p.scales(), p.dnorm()), Andreas Ruckstuhl [ctb] (Functions: p.arrows(), p.profileTraces(), p.res.2x()), Christian Keller [ctb] (Functions: histBxp(), p.tachoPlot()), Kjetil Halvorsen [ctb] (Functions: KSd(), ecdf.ksCI()), Alain Hauser [ctb] (Functions: cairoSwd(), is.whole(), toLatex.numeric()*), Christoph Buser [ctb] (to function Duplicated()), Lorenz Gygax [ctb] (to function p.res.2fact()), Bill Venables [ctb] (Functions: empty.dimnames(), primes()), Tony Plate [ctb] (to inv.seq()), Isabelle Fl<fc>ckiger [ctb], Marcel Wolbers [ctb], Markus Keller [ctb], Sandrine Dudoit [ctb], Jane Fridlyand [ctb], Greg Snow [ctb] (to loessDemo()), Henrik Aa. Nielsen [ctb] (to loessDemo()), Vincent Carey [ctb], Ben Bolker [ctb], Philippe Grosjean [ctb], Fr<e9>d<e9>ric Ibanez [ctb], Caterina Savi [ctb], Charles Geyer [ctb], Jens Oehlschl<e4>gel [ctb]
Initial release
2021-04-03

We don't support your browser anymore

Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.