Semi-Primes
September 6, 2019
If you have a function that lists primes, it is easy to make a function that lists semi-primes:
(define (semi-primes n)
(let ((ps (primes (quotient n 2))))
(let loop ((ps ps) (qs (cdr ps)) (ss (list)))
(cond ((< n (* (car ps) (cadr ps))) (sort < ss))
((or (null? qs) (< n (* (car ps) (car qs))))
(loop (cdr ps) (cddr ps) ss))
(else (loop ps (cdr qs) (cons (* (car ps) (car qs)) ss)))))))
The program iterates over the primes ps, looping through the remaining primes at each step to form semi-primes, stopping when the product exceeds the desired n. Here we calculate the number of semi-primes at each power of ten (A036351).
> (do ((k 1 (+ k 1))) ((= k 8))
(display k) (display #\tab)
(display (length (semi-primes (expt 10 k))))
(newline))
1 2
2 30
3 288
4 2600
5 23313
6 209867
7 1903878
You can run the program at https://ideone.com/jH2sF6.
In Python. The semiprimes are generated in order.
from heapq import merge from itertools import takewhile def semiprimes(limit): prs = list(primes(limit // 2)) seqs = [list(takewhile(lambda p: p <= limit, (a*b for b in prs[i+1:]))) for i, a in enumerate(prs)] return list(merge(*seqs)) print(semiprimes(25)) # -> [6, 10, 14, 15, 21, 22] print(len(semiprimes(10**4))) # -> 2600An faster version in Python. The Python sort function is very efficient for a sequence of runs of already sorted items, so there is no point in using merge.
def semiprimes4(limit): prs = list(primes(limit // 2)) semis = [] for i, a in enumerate(prs): if a * prs[i+1] > limit: # saves a lot of time break for b in prs[i+1:]: ab = a * b if ab > limit: break semis.append(ab) return sorted(semis)bool isPrime(int n)
{
if (n <= 1)
return false;
Here’s some Haskell, making use of lazy lists:
All those merges build up in memory, but we can find eg. 5111663, the 1000000th semiprime without too much trouble in 30 seconds or so.
If we include the squares of primes in the semiprimes, the code is even simpler: