17 ms·
My favorite prime number generator
- anaphor 3y agoThere's a functional / lazy version I like that uses a priority queue in order to get the next prime number https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf It basically uses the priority queue to store all the composites up to a certain point and then you can assume the next one is prime. I had some (admittedly not well written) code for it too https://gist.github.com/weskerfoot/4699275 https://gist.github.com/weskerfoot/4699275 It can also incorporates the wheel factorization optimization mentioned in this article.
- oxxoxoxooo 3y agoAnd another one: The Genuine Sieve of Eratosthenes https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf
- tzot 3y agoThe fastest variant I know of this algorithm is the one found at https://stackoverflow.com/a/3796442/6899 https://stackoverflow.com/a/3796442/6899 : import itertools as it def erat3( ): D = { 9: 3, 25: 5 } yield 2 yield 3 yield 5 MASK= 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0, MODULOS= frozenset( (1, 7, 11, 13, 17, 19, 23, 29) ) for q in it.compress( it.islice(it.count(7), 0, None, 2), it.cycle(MASK)): p = D.pop(q, None) if p is None: D[q*q] = q yield q else: x = q + 2*p while x in D or (x%30) not in MODULOS: x += 2*p D[x] = p
- svat 3y agoThis adds two things over the OP's initial version, both mentioned in the post: wheel factorization (with 30), and the "linear probing" to avoid lists. It's possible to keep improving (up to a point) by using e.g. wheel with 210 instead of 30, and so on, if one doesn't mind hard-coding more and more.
- simlevesque 3y agoisn't hardcoding 2 3 and 5 a bit of a cheat ? It seems like you could hardcode it to yield the X next primes and it'd be faster.
- stouset 3y agoThe point isn’t hardcoding those numbers. The point is that you don’t iterate on any multiples of those numbers (which we know are composite) and doing allows us to eliminate 2/3 of the iterations. It’s an extension of enumerating over only odd numbers, which eliminates 1/2 of the iterations. Hardcoding those values is a semi-requirement of the above technique.
- simlevesque 3y agoThank you, it makes more sense now.
- Exuma 3y agoHe said he has "used this code many times"... have you ever needed to find prime numbers in your daily programming?
- thrownblown 3y ago[flagged]
- thrownblown 3y agoto those who flagged my "do you even crypto bro?" comment: even a little research into cryptography will inform you that prime numbers are used heavily in encryption algorithms. https://stackoverflow.com/questions/439870/why-are-primes-important-in-cryptography https://stackoverflow.com/questions/439870/why-are-primes-im...
- Exuma 3y agoDon't try to pull a "even a little research" bullshit comment on me. "even a little research" will tell you people shouldn't roll your own cryptography and to just use a library. So, this still doesn't answer my question.
- thrownblown 3y agoyou'll never be in cicada 3301 with that kinda attitude my guy! but to counter your point: people do have to write those libraries that you depend on. others have to validate them, and yet still others are doing their best efforts to invalidate them and attack them. but you are correct that most people should not roll their own cryptography suites.
- krnaveen14 3y agoReminds me of the time where I've implemented highly optimised version of Sieve of Eratosthenes in Java, just for a Code Golf Challenge. Memory Usage is greatly reduced by using BitSet. https://github.com/krnaveen14/PrimeNumbers/blob/master/src/main/java/com/procreators/primenumber/CustomSieveOfEratosthenes.java https://github.com/krnaveen14/PrimeNumbers/blob/master/src/m... https://codegolf.stackexchange.com/a/66915/41984 https://codegolf.stackexchange.com/a/66915/41984 OMG, it's been 8 years already!
- omoikane 3y agoMy favorite prime number generator is the undocumented __next_prime(): https://github.com/llvm-mirror/libcxx/blob/78d6a7767ed57b50122a161b91f59f19c9bd0d19/src/hash.cpp#L170 https://github.com/llvm-mirror/libcxx/blob/78d6a7767ed57b501... There is no good reason to use this one except in a code golf environment that includes all headers by default, which is where I learned about it.
- manicennui 3y agoI've always enjoyed the Perl regex solution.
- SeanLuke 3y agoBelieve it or not, there exists a closed-form equation which computes the nth prime number. It's Willans's formula: p_n = 1 + Sum_{i=1}^{2^n} Floor((n/(Sum_{j=1}^{i} Floor((Cos((j-1)!+1)/j) Pi)^2)))^(1/n)) https://en.wikipedia.org/wiki/Formula_for_primes https://en.wikipedia.org/wiki/Formula_for_primes https://www.youtube.com/watch?v=j5s0h42GfvM https://www.youtube.com/watch?v=j5s0h42GfvM It's not fast but it really does work.
- svat 3y agoThat video is actually great! The last time this formula came up, I added some references to the Wikipedia article (current version: https://en.wikipedia.org/w/index.php?title=Formula_for_primes&oldid=1170749145#Formulas_based_on_Wilson's_theorem) https://en.wikipedia.org/w/index.php?title=Formula_for_prime... about the "worthlessness" of such formulas, but it doesn't deter people from being interested in that formula… as long as they don't imagine it's useful in any way, it's fine to enjoy its perverse cleverness, I suppose. The video does well to quote Wilf on “the disease was preferable to the cure” and to end with “Could it [C. P. Willans] be a pseudonym for someone who didn't want to tarnish their reputation by writing about useless prime-generating formulas?” :-) Dudley (who mentions a large list of other such formulas) also has a nice conclusion: > The conclusion to be drawn from all this, I think, is that formulas for formulas' sake do not advance the mathematical enterprise. Formulas should be useful. If not, they should be astounding, elegant, enlightening, simple, or have some other redeeming value.
- theshrike79 3y agoI just use the Prime Number Shitting Bear personally: http://alpha61.com/primenumbershittingbear/ http://alpha61.com/primenumbershittingbear/
- bjoli 3y agoThis is a paper which describes a fast lazy sieve of ersthostenes: https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf https://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf
- intc 3y agoSimple C based approach on Sieve of Eratosthenes (by yours truly): * https://github.com/intc/c-notes/blob/master/threads/4-thrd_eratosthenes/erat.c https://github.com/intc/c-notes/blob/master/threads/4-thrd_e...
- jacknews 3y agoNice, I like the optimization replacing the list, and just hunting for the next 'free' composite for the prime. Teaching my kids coding, we also arrived at most of the other optimizations, in Janet: (defn primes-sieve [] (fiber/new (fn [] (yield 2) # the only even prime, yield it explicitly. (let [sieve (table/new (math/pow 2 12))] (loop [n :range [3 nil 2]] # only odd numbers (if-let [current-composite (sieve n)] (do (each prime-factor current-composite (let [next-composite-key (+ n (* 2 prime-factor))] # add 2*prime-factor, else composite would be even (if-let [next-composite (sieve next-composite-key)] (array/push next-composite prime-factor) (put sieve next-composite-key @[prime-factor])))) (put sieve n nil)) (do (yield n) (put sieve (* n n) @[n])) # next composite is n^2. Anything less is composite a smaller prime, already in the sieve ))))))
- tromp 3y agoMy favorite prime number generator is ┬─┬───────────────────────────────────┬──── ────┬ │ │ ┬─┬ ────┬─┬────────────────────── ┼───┬ ┬───┼ │ │ └─┤ ────┼─┼───────────────────┬── ┼───┼ │ ┬ │ │ │ │ ┬───┼─┼───────────────────┼── ┼─┬─┼ │ ┼ │ │ │ │ │ ─ ┼─┼─┬─────┬──── ──────┼─┬ │ ├─┘ └─┤ │ │ │ │ │ ┬ └─┤ │ ┬─┬ ┼─┬─┬ ──┬───┼─┼ ├─┘ ├─┘ │ │ │ └─┤ └─┤ └─┤ │ ├─┘ ──┼─┬─┼─┼ │ │ │ │ │ │ │ │ ├─┘ ┬─┼─┼─┼─┼ │ │ │ │ │ │ │ ├─┘ └─┤ │ ├─┘ │ │ │ │ │ │ └───┤ │ ├─┘ │ │ │ │ │ │ │ ├─┘ │ │ │ │ │ │ ├─────────┘ │ │ │ │ │ ├─────────┘ │ │ │ │ └───┤ │ │ │ │ ├───────────────────────────┘ │ │ ├───────┘ │ └─┤ │ └───────────────────────────────────────────┘ which graphically depicts the 167-bit lambda term λ 1 (1 ((λ 1 1) (λ λ λ 1 (λ λ 1) ((λ 4 4 1 ((λ 1 1) (λ 2 (1 1)))) (λ λ λ λ 1 3 (2 (6 4))))) (λ λ λ 4 (1 3)))) (λ λ 1 (λ λ 2) 2) which generates the characteristic sequence of primes [1] [2]. [1] https://tromp.github.io/cl/cl.html https://tromp.github.io/cl/cl.html [2] https://github.com/tromp/AIT/blob/master/characteristic_sequences/primes.lam https://github.com/tromp/AIT/blob/master/characteristic_sequ...
- hgsgm 3y agoWhat's the "characteristic sequence"?
- hidroto 3y agoits a list with each prime index contains a one/true value
- cobaltoxide 3y agoCould you explain it? Those two links are Greek to me.
- tromp 3y agoThe first link starts out as: > On Pi day 2023 I gave an online talk about AIT and BLC based on these slides. Perhaps watching the linked-to talk will somewhat clarify matters.
- naillo 3y agoThis person has an incredible work ethic
- shaftoe444 3y agoMy favourite prime number generator is this one using go channels to create filters. Not efficient but fun. https://go.dev/play/p/9U22NfrXeq https://go.dev/play/p/9U22NfrXeq
- tantalor 3y agoMy favorite prime number generator is http://alpha61.com/primenumbershittingbear/ http://alpha61.com/primenumbershittingbear/
- proactivesvcs 3y agoIt doesn't need Flash Player any more? It's the year 2023 and I can run the Prime Number Shitting Bear again? Hallelujah!
- tyingq 3y agoIt is still flash. Appears they incorporated aflax, a "run most of flash in js" library. https://github.com/aptana/aflax https://github.com/aptana/aflax
- deleted 3y ago[deleted]
- max_ 3y agoIn J generating prime numbers is simply p:21, that gives you 79 which is the 21st prime number. It is so fast, I wonder what the algorithm behind it is. [0]: https://code.jsoftware.com/wiki/Vocabulary/pco https://code.jsoftware.com/wiki/Vocabulary/pco
- jonahx 3y ago> I wonder what the algorithm behind it is. Noted in the details section of the link you provided: "Primality testing on numbers larger than (2^31) uses the probabilistic Miller-Rabin algorithm."
- noman-land 3y agoHow high can you go with this?
- max_ 3y agoI have J on my android phone, and the highest I could go was p: 100999999, the prime numbers is 2059519669
- nemo8551 3y agoCan I just say how surprised I am the many people have a favourite prime number generator. I can’t say I have a favourite but I do enjoy multiplying two random prime numbers with each other and using the random result when I need a random number with added random.
- cobaltoxide 3y agoThat doesn't seem like a good way to generate "random" numbers, since it will only produce numbers with exactly two prime factors.
- amenghra 3y agoCorrect. It will almost never produce an even number, so there's already one bit of entropy lost. It will also never produce lots of other numbers.
- nemo8551 3y agoThanks all, I hadn’t considered that but we do only use it for silly things and nothing work related but it’s worth a consideration.
- andrelaszlo 3y agoA way to add even more random is to generate a list of N random primes, then shuffling it M times before multiplying the elements. N and M are of course random positive integers. To add maximum random, sometimes I like to return a string with a random Gertrude Stein quote, segfault, or make the PC speaker beep instead of returning a number. So random!
- ChrisCinelli 3y agoThe Sieve of Atkin is VERY fast at generating primes between 1 and n. It asymptotically speeds up the process of generating prime numbers. It has time complexity O(n/(log log n)) https://www.baeldung.com/cs/prime-number-algorithms https://www.baeldung.com/cs/prime-number-algorithms
- anonymoushn 3y agoYour link seems to say it has time complexity O(n)?
- hgsgm 3y agoIf you count everything carefully, it's O(n / log log n) https://www.ams.org/journals/mcom/2004-73-246/S0025-5718-03-01501-1/S0025-5718-03-01501-1.pdf https://www.ams.org/journals/mcom/2004-73-246/S0025-5718-03-... They are essentially the same in all conceivable practical cases in this Universe. log log (atoms in universe) < 100
- anonymoushn 3y agoThis seems to be a matter of the above source reporting a variation of the algorithm that is asymptotically slower than the algorithm given in the paper *and* uses asymptotically more memory. Thanks for posting the paper!
- svat 3y agoThis is really elegant! To make the explanation more concrete, here's a snapshot of the algorithm's state (the dict D) after it has processed numbers up to 10 (inclusive): 12: [3, 2] 25: [5] 49: [7] So each prime up to 10 is present exactly once. When we process 12, we move 2 to the list for 14, and 3 to the list for 15. The Sorenson paper linked in the post is also a great read, I remember reading it a couple of years ago and it is very clear. It adds many further optimizations on top of standard sieve of Eratosthenes (e.g. wheel factorization is a generalization of the idea of considering only numbers in {1, 5} mod 6, or {1, 7, 11, 13, 17, 19, 23, 29} mod 30, or …), and segmenting is what you need to do to avoid running out of memory. The fastest such generator (using sieve-based methods) seems to be https://github.com/kimwalisch/primesieve https://github.com/kimwalisch/primesieve — it implements all those optimizations and more. However, if you just want the nth prime, it's even faster to compute π(n) by combinatorial/analytic methods, instead of a sieve: the same author has https://github.com/kimwalisch/primecount https://github.com/kimwalisch/primecount for that. (Try the ten-trillionth prime number.)
- scythe 3y agoThe memory usage can in fact be optimized further with a small change: https://gist.github.com/scythe/5fb962722934c58c60430180beab86ed https://gist.github.com/scythe/5fb962722934c58c60430180beab8... In the spirit of the blog post, I'll let you guess how it works :p
- svat 3y agoWow this is a great optimization! You're under-selling it :) It changes the memory requirement from O(π(n)) to O(π(√n)), e.g. after processing numbers up to 100, the internal state in the OP's algorithm is: q = 101 D = {'102': [3, 2], '105': [7, 5], '121': [11], '169': [13], '289': [17], '361': [19], '529': [23], '841': [29], '961': [31], '1369': [37], '1681': [41], '1849': [43], '2209': [47], '2809': [53], '3481': [59], '3721': [61], '4489': [67], '5041': [71], '5329': [73], '6241': [79], '6889': [83], '7921': [89], '9409': [97]} while with your "small change", the state is just: q = 101 S = 121 r = 11 D = {'102': [3, 2], '105': [7, 5]}
- philipLutz 3y agoI enjoyed this, but I am a little stuck on how stuck the author is on yield. If one ends up creating a large list of prime numbers with this generator, why not just return a large list of prime numbers? The author clearly notes how optimizing Python code is "odd" when you can more easily switch languages, but languages like C/C++/Rust would not have pretty code because it does not have yield. It is obvious that returning an array of prime numbers is what one would do in C/C++/Rust. I understand that not using yield would mean it is not strictly a "generator", but the spirit of this problem is generating the sequence of primes (unless I am missing something).
- jameshart 3y agoIf the problem is outputting the primes, then storing them in a list is not a necessary part of that problem. When writing ‘fizz buzz’, would you expect to return an array of 100 strings, some of which are numbers, some ‘fizz’es and some ‘buzz’es?
- vlovich123 3y ago> but languages like C/C++/Rust would not have pretty code because it does not have yield Maybe I'm misunderstanding you, but Rust has native generators with yield: https://doc.rust-lang.org/beta/unstable-book/language-features/generators.html https://doc.rust-lang.org/beta/unstable-book/language-featur... C++ also does using coroutines + co_yield (https://en.cppreference.com/w/cpp/language/coroutines#co_yield https://en.cppreference.com/w/cpp/language/coroutines#co_yie...) or iterators (old syntax).
- masklinn 3y ago> If one ends up creating a large list of prime numbers with this generator But one does not always end up creating a large list of prime numbers, do they? > why not just return a large list of prime numbers? Because that requires upfront knowledge of the filtering factor or absence thereof. A generator means you don't care, you generate an infinite sequence and the consumer is free to do whatever they want. > It is obvious that returning an array of prime numbers is what one would do in C/C++/Rust. It's certainly not obvious for Rust.
- redocneknurd 3y agoWhat are you using your prime numbers for?
- earleybird 3y agoPaaS - Primes as a Service; spun out of the IaaS, Integers as a Service
- dvh 3y ago{(~T∊T∘.×T)/T←1↓⍳⍵}
- jodrellblank 3y ago{(~T∊T∘.×T)/T←1↓⍳⍵} 25 {T~∘∊∘.×⍨T←1↓⍳⍵} 25 ((⊢~∘∊∘.×⍨)1↓⍳) 25
- ducttapecrown 3y agoWhat language is that? AWK??
- paulddraper 3y agoIt's got funny symbols, so I assume APL. [1] [1] https://en.wikipedia.org/wiki/APL_(programming_language) https://en.wikipedia.org/wiki/APL_(programming_language)
- jodrellblank 3y agoIt's an APL sieve, for an input of 15 (primes up to 15) it does these steps: First 15 integers: ⍳ 15 ┌→──────────────────────────────────┐ │1 2 3 4 5 6 7 8 9 10 11 12 13 14 15│ └~──────────────────────────────────┘ 1 element dropped from the front of that: 1↓⍳15 ┌→────────────────────────────────┐ │2 3 4 5 6 7 8 9 10 11 12 13 14 15│ └~────────────────────────────────┘ That array stored in variable T: T←1↓⍳15 A table of outer-product multiplication, T by T: T∘.×T ┌→─────────────────────────────────────────────────┐ ↓ 4 6 8 10 12 14 16 18 20 22 24 26 28 30│ │ 6 9 12 15 18 21 24 27 30 33 36 39 42 45│ │ 8 12 16 20 24 28 32 36 40 44 48 52 56 60│ │10 15 20 25 30 35 40 45 50 55 60 65 70 75│ │12 18 24 30 36 42 48 54 60 66 72 78 84 90│ │14 21 28 35 42 49 56 63 70 77 84 91 98 105│ │16 24 32 40 48 56 64 72 80 88 96 104 112 120│ │18 27 36 45 54 63 72 81 90 99 108 117 126 135│ │20 30 40 50 60 70 80 90 100 110 120 130 140 150│ │22 33 44 55 66 77 88 99 110 121 132 143 154 165│ │24 36 48 60 72 84 96 108 120 132 144 156 168 180│ │26 39 52 65 78 91 104 117 130 143 156 169 182 195│ │28 42 56 70 84 98 112 126 140 154 168 182 196 210│ │30 45 60 75 90 105 120 135 150 165 180 195 210 225│ └~─────────────────────────────────────────────────┘ The elements of T which exist (are found) in the multiplication table, boolean result: T∊T∘.×T ┌→──────────────────────────┐ │0 0 1 0 1 0 1 1 1 0 1 0 1 1│ └~──────────────────────────┘ Showing that together with the input, 4 6 8 are found, 2 3 5 are not: T,[0.5]T∊T∘.×T ┌→────────────────────────────────┐ ↓2 3 4 5 6 7 8 9 10 11 12 13 14 15│ │0 0 1 0 1 0 1 1 1 0 1 0 1 1│ └~────────────────────────────────┘ Numbers not found in the multiplication table are prime, so invert the test result so primes get 1 and composites get 0: ~T∊T∘.×T ┌→──────────────────────────┐ │1 1 0 1 0 1 0 0 0 1 0 1 0 0│ └~──────────────────────────┘ Use that to 'compress' the input to keep things where the 1 positions are and drop things where the 0s are: (~T∊T∘.×T)/T ┌→────────────┐ │2 3 5 7 11 13│ └~────────────┘ You can play with it at https://tryapl.org/ https://tryapl.org/ using the bar at the top to enter the symbols.
- tzs 3y agoMy favorite is this FRACTRAN program: (17/91, 78/85, 19/51, 23/38, 29/33, 77/29, 1/17, 11/13, 13/11, 15/2, 1/7, 55/1) with input 2. It doesn't quite produce primes directly. It gives a sequence of integers, some of which are powers of 2. Those powers of 2 are 2^2, 2^3, 2^5, 2^7, 2^11, ..., i.e., the sequence of 2 to the power of primes. For those who have never seen FRACTRAN it is an esoteric programming language invented by John Conway. It is Turing-complete. Here's how you execute a FRACTRAN program, with input N. 1. Step through the list of fractions in order until you find a fraction f such that fN is an integer or you run out of fractions. 2. If you run out of fractions halt. 3. Output fN. 4. Replace N with fN. 5. Goto step 1. Wikipedia gives some sample programs and explains in detail how the heck FRACTRAN can compute: https://en.wikipedia.org/wiki/FRACTRAN https://en.wikipedia.org/wiki/FRACTRAN
- rullelito 3y agoMy favorite one: d = defaultdict(list) # Magic for i in range(2, 10000): for y in i in d and d.pop(i) or [i]: d[i+y] += (y,) # Not sure how to print this in a succinct way print(sorted(list(set().union(*d.values()))))
- paulddraper 3y agoBeautiful from collections import defaultdict d = defaultdict(list) for i in range(2, 10000): for y in d.pop(i, [i]): d[i + y].append(y) for y in sorted(y for x in d.values() for y in x): print(y)