109c434b8SThomas Gleixner // SPDX-License-Identifier: GPL-2.0-only
2*9ac17575SChristophe JAILLET #define pr_fmt(fmt) "prime numbers: " fmt
32c64e9cbSAndy Shevchenko
42c64e9cbSAndy Shevchenko #include <linux/module.h>
52c64e9cbSAndy Shevchenko #include <linux/mutex.h>
62c64e9cbSAndy Shevchenko #include <linux/prime_numbers.h>
72c64e9cbSAndy Shevchenko #include <linux/slab.h>
82c64e9cbSAndy Shevchenko
92c64e9cbSAndy Shevchenko struct primes {
102c64e9cbSAndy Shevchenko struct rcu_head rcu;
112c64e9cbSAndy Shevchenko unsigned long last, sz;
122c64e9cbSAndy Shevchenko unsigned long primes[];
132c64e9cbSAndy Shevchenko };
142c64e9cbSAndy Shevchenko
152c64e9cbSAndy Shevchenko #if BITS_PER_LONG == 64
162c64e9cbSAndy Shevchenko static const struct primes small_primes = {
172c64e9cbSAndy Shevchenko .last = 61,
182c64e9cbSAndy Shevchenko .sz = 64,
192c64e9cbSAndy Shevchenko .primes = {
202c64e9cbSAndy Shevchenko BIT(2) |
212c64e9cbSAndy Shevchenko BIT(3) |
222c64e9cbSAndy Shevchenko BIT(5) |
232c64e9cbSAndy Shevchenko BIT(7) |
242c64e9cbSAndy Shevchenko BIT(11) |
252c64e9cbSAndy Shevchenko BIT(13) |
262c64e9cbSAndy Shevchenko BIT(17) |
272c64e9cbSAndy Shevchenko BIT(19) |
282c64e9cbSAndy Shevchenko BIT(23) |
292c64e9cbSAndy Shevchenko BIT(29) |
302c64e9cbSAndy Shevchenko BIT(31) |
312c64e9cbSAndy Shevchenko BIT(37) |
322c64e9cbSAndy Shevchenko BIT(41) |
332c64e9cbSAndy Shevchenko BIT(43) |
342c64e9cbSAndy Shevchenko BIT(47) |
352c64e9cbSAndy Shevchenko BIT(53) |
362c64e9cbSAndy Shevchenko BIT(59) |
372c64e9cbSAndy Shevchenko BIT(61)
382c64e9cbSAndy Shevchenko }
392c64e9cbSAndy Shevchenko };
402c64e9cbSAndy Shevchenko #elif BITS_PER_LONG == 32
412c64e9cbSAndy Shevchenko static const struct primes small_primes = {
422c64e9cbSAndy Shevchenko .last = 31,
432c64e9cbSAndy Shevchenko .sz = 32,
442c64e9cbSAndy Shevchenko .primes = {
452c64e9cbSAndy Shevchenko BIT(2) |
462c64e9cbSAndy Shevchenko BIT(3) |
472c64e9cbSAndy Shevchenko BIT(5) |
482c64e9cbSAndy Shevchenko BIT(7) |
492c64e9cbSAndy Shevchenko BIT(11) |
502c64e9cbSAndy Shevchenko BIT(13) |
512c64e9cbSAndy Shevchenko BIT(17) |
522c64e9cbSAndy Shevchenko BIT(19) |
532c64e9cbSAndy Shevchenko BIT(23) |
542c64e9cbSAndy Shevchenko BIT(29) |
552c64e9cbSAndy Shevchenko BIT(31)
562c64e9cbSAndy Shevchenko }
572c64e9cbSAndy Shevchenko };
582c64e9cbSAndy Shevchenko #else
592c64e9cbSAndy Shevchenko #error "unhandled BITS_PER_LONG"
602c64e9cbSAndy Shevchenko #endif
612c64e9cbSAndy Shevchenko
622c64e9cbSAndy Shevchenko static DEFINE_MUTEX(lock);
632c64e9cbSAndy Shevchenko static const struct primes __rcu *primes = RCU_INITIALIZER(&small_primes);
642c64e9cbSAndy Shevchenko
652c64e9cbSAndy Shevchenko static unsigned long selftest_max;
662c64e9cbSAndy Shevchenko
slow_is_prime_number(unsigned long x)672c64e9cbSAndy Shevchenko static bool slow_is_prime_number(unsigned long x)
682c64e9cbSAndy Shevchenko {
692c64e9cbSAndy Shevchenko unsigned long y = int_sqrt(x);
702c64e9cbSAndy Shevchenko
712c64e9cbSAndy Shevchenko while (y > 1) {
722c64e9cbSAndy Shevchenko if ((x % y) == 0)
732c64e9cbSAndy Shevchenko break;
742c64e9cbSAndy Shevchenko y--;
752c64e9cbSAndy Shevchenko }
762c64e9cbSAndy Shevchenko
772c64e9cbSAndy Shevchenko return y == 1;
782c64e9cbSAndy Shevchenko }
792c64e9cbSAndy Shevchenko
slow_next_prime_number(unsigned long x)802c64e9cbSAndy Shevchenko static unsigned long slow_next_prime_number(unsigned long x)
812c64e9cbSAndy Shevchenko {
822c64e9cbSAndy Shevchenko while (x < ULONG_MAX && !slow_is_prime_number(++x))
832c64e9cbSAndy Shevchenko ;
842c64e9cbSAndy Shevchenko
852c64e9cbSAndy Shevchenko return x;
862c64e9cbSAndy Shevchenko }
872c64e9cbSAndy Shevchenko
clear_multiples(unsigned long x,unsigned long * p,unsigned long start,unsigned long end)882c64e9cbSAndy Shevchenko static unsigned long clear_multiples(unsigned long x,
892c64e9cbSAndy Shevchenko unsigned long *p,
902c64e9cbSAndy Shevchenko unsigned long start,
912c64e9cbSAndy Shevchenko unsigned long end)
922c64e9cbSAndy Shevchenko {
932c64e9cbSAndy Shevchenko unsigned long m;
942c64e9cbSAndy Shevchenko
952c64e9cbSAndy Shevchenko m = 2 * x;
962c64e9cbSAndy Shevchenko if (m < start)
972c64e9cbSAndy Shevchenko m = roundup(start, x);
982c64e9cbSAndy Shevchenko
992c64e9cbSAndy Shevchenko while (m < end) {
1002c64e9cbSAndy Shevchenko __clear_bit(m, p);
1012c64e9cbSAndy Shevchenko m += x;
1022c64e9cbSAndy Shevchenko }
1032c64e9cbSAndy Shevchenko
1042c64e9cbSAndy Shevchenko return x;
1052c64e9cbSAndy Shevchenko }
1062c64e9cbSAndy Shevchenko
expand_to_next_prime(unsigned long x)1072c64e9cbSAndy Shevchenko static bool expand_to_next_prime(unsigned long x)
1082c64e9cbSAndy Shevchenko {
1092c64e9cbSAndy Shevchenko const struct primes *p;
1102c64e9cbSAndy Shevchenko struct primes *new;
1112c64e9cbSAndy Shevchenko unsigned long sz, y;
1122c64e9cbSAndy Shevchenko
1132c64e9cbSAndy Shevchenko /* Betrand's Postulate (or Chebyshev's theorem) states that if n > 3,
1142c64e9cbSAndy Shevchenko * there is always at least one prime p between n and 2n - 2.
1152c64e9cbSAndy Shevchenko * Equivalently, if n > 1, then there is always at least one prime p
1162c64e9cbSAndy Shevchenko * such that n < p < 2n.
1172c64e9cbSAndy Shevchenko *
1182c64e9cbSAndy Shevchenko * http://mathworld.wolfram.com/BertrandsPostulate.html
1192c64e9cbSAndy Shevchenko * https://en.wikipedia.org/wiki/Bertrand's_postulate
1202c64e9cbSAndy Shevchenko */
1212c64e9cbSAndy Shevchenko sz = 2 * x;
1222c64e9cbSAndy Shevchenko if (sz < x)
1232c64e9cbSAndy Shevchenko return false;
1242c64e9cbSAndy Shevchenko
1252c64e9cbSAndy Shevchenko sz = round_up(sz, BITS_PER_LONG);
1262c64e9cbSAndy Shevchenko new = kmalloc(sizeof(*new) + bitmap_size(sz),
1272c64e9cbSAndy Shevchenko GFP_KERNEL | __GFP_NOWARN);
1282c64e9cbSAndy Shevchenko if (!new)
1292c64e9cbSAndy Shevchenko return false;
1302c64e9cbSAndy Shevchenko
1312c64e9cbSAndy Shevchenko mutex_lock(&lock);
1322c64e9cbSAndy Shevchenko p = rcu_dereference_protected(primes, lockdep_is_held(&lock));
1332c64e9cbSAndy Shevchenko if (x < p->last) {
1342c64e9cbSAndy Shevchenko kfree(new);
1352c64e9cbSAndy Shevchenko goto unlock;
1362c64e9cbSAndy Shevchenko }
1372c64e9cbSAndy Shevchenko
1382c64e9cbSAndy Shevchenko /* Where memory permits, track the primes using the
1392c64e9cbSAndy Shevchenko * Sieve of Eratosthenes. The sieve is to remove all multiples of known
1402c64e9cbSAndy Shevchenko * primes from the set, what remains in the set is therefore prime.
1412c64e9cbSAndy Shevchenko */
1422c64e9cbSAndy Shevchenko bitmap_fill(new->primes, sz);
1432c64e9cbSAndy Shevchenko bitmap_copy(new->primes, p->primes, p->sz);
1442c64e9cbSAndy Shevchenko for (y = 2UL; y < sz; y = find_next_bit(new->primes, sz, y + 1))
1452c64e9cbSAndy Shevchenko new->last = clear_multiples(y, new->primes, p->sz, sz);
1462c64e9cbSAndy Shevchenko new->sz = sz;
1472c64e9cbSAndy Shevchenko
1482c64e9cbSAndy Shevchenko BUG_ON(new->last <= x);
1492c64e9cbSAndy Shevchenko
1502c64e9cbSAndy Shevchenko rcu_assign_pointer(primes, new);
1512c64e9cbSAndy Shevchenko if (p != &small_primes)
1522c64e9cbSAndy Shevchenko kfree_rcu((struct primes *)p, rcu);
1532c64e9cbSAndy Shevchenko
1542c64e9cbSAndy Shevchenko unlock:
1552c64e9cbSAndy Shevchenko mutex_unlock(&lock);
1562c64e9cbSAndy Shevchenko return true;
1572c64e9cbSAndy Shevchenko }
1582c64e9cbSAndy Shevchenko
free_primes(void)1592c64e9cbSAndy Shevchenko static void free_primes(void)
1602c64e9cbSAndy Shevchenko {
1612c64e9cbSAndy Shevchenko const struct primes *p;
1622c64e9cbSAndy Shevchenko
1632c64e9cbSAndy Shevchenko mutex_lock(&lock);
1642c64e9cbSAndy Shevchenko p = rcu_dereference_protected(primes, lockdep_is_held(&lock));
1652c64e9cbSAndy Shevchenko if (p != &small_primes) {
1662c64e9cbSAndy Shevchenko rcu_assign_pointer(primes, &small_primes);
1672c64e9cbSAndy Shevchenko kfree_rcu((struct primes *)p, rcu);
1682c64e9cbSAndy Shevchenko }
1692c64e9cbSAndy Shevchenko mutex_unlock(&lock);
1702c64e9cbSAndy Shevchenko }
1712c64e9cbSAndy Shevchenko
1722c64e9cbSAndy Shevchenko /**
1732c64e9cbSAndy Shevchenko * next_prime_number - return the next prime number
1742c64e9cbSAndy Shevchenko * @x: the starting point for searching to test
1752c64e9cbSAndy Shevchenko *
1762c64e9cbSAndy Shevchenko * A prime number is an integer greater than 1 that is only divisible by
1772c64e9cbSAndy Shevchenko * itself and 1. The set of prime numbers is computed using the Sieve of
1782c64e9cbSAndy Shevchenko * Eratoshenes (on finding a prime, all multiples of that prime are removed
1792c64e9cbSAndy Shevchenko * from the set) enabling a fast lookup of the next prime number larger than
1802c64e9cbSAndy Shevchenko * @x. If the sieve fails (memory limitation), the search falls back to using
1812c64e9cbSAndy Shevchenko * slow trial-divison, up to the value of ULONG_MAX (which is reported as the
1822c64e9cbSAndy Shevchenko * final prime as a sentinel).
1832c64e9cbSAndy Shevchenko *
1842c64e9cbSAndy Shevchenko * Returns: the next prime number larger than @x
1852c64e9cbSAndy Shevchenko */
next_prime_number(unsigned long x)1862c64e9cbSAndy Shevchenko unsigned long next_prime_number(unsigned long x)
1872c64e9cbSAndy Shevchenko {
1882c64e9cbSAndy Shevchenko const struct primes *p;
1892c64e9cbSAndy Shevchenko
1902c64e9cbSAndy Shevchenko rcu_read_lock();
1912c64e9cbSAndy Shevchenko p = rcu_dereference(primes);
1922c64e9cbSAndy Shevchenko while (x >= p->last) {
1932c64e9cbSAndy Shevchenko rcu_read_unlock();
1942c64e9cbSAndy Shevchenko
1952c64e9cbSAndy Shevchenko if (!expand_to_next_prime(x))
1962c64e9cbSAndy Shevchenko return slow_next_prime_number(x);
1972c64e9cbSAndy Shevchenko
1982c64e9cbSAndy Shevchenko rcu_read_lock();
1992c64e9cbSAndy Shevchenko p = rcu_dereference(primes);
2002c64e9cbSAndy Shevchenko }
2012c64e9cbSAndy Shevchenko x = find_next_bit(p->primes, p->last, x + 1);
2022c64e9cbSAndy Shevchenko rcu_read_unlock();
2032c64e9cbSAndy Shevchenko
2042c64e9cbSAndy Shevchenko return x;
2052c64e9cbSAndy Shevchenko }
2062c64e9cbSAndy Shevchenko EXPORT_SYMBOL(next_prime_number);
2072c64e9cbSAndy Shevchenko
2082c64e9cbSAndy Shevchenko /**
2092c64e9cbSAndy Shevchenko * is_prime_number - test whether the given number is prime
2102c64e9cbSAndy Shevchenko * @x: the number to test
2112c64e9cbSAndy Shevchenko *
2122c64e9cbSAndy Shevchenko * A prime number is an integer greater than 1 that is only divisible by
2132c64e9cbSAndy Shevchenko * itself and 1. Internally a cache of prime numbers is kept (to speed up
2142c64e9cbSAndy Shevchenko * searching for sequential primes, see next_prime_number()), but if the number
2152c64e9cbSAndy Shevchenko * falls outside of that cache, its primality is tested using trial-divison.
2162c64e9cbSAndy Shevchenko *
2172c64e9cbSAndy Shevchenko * Returns: true if @x is prime, false for composite numbers.
2182c64e9cbSAndy Shevchenko */
is_prime_number(unsigned long x)2192c64e9cbSAndy Shevchenko bool is_prime_number(unsigned long x)
2202c64e9cbSAndy Shevchenko {
2212c64e9cbSAndy Shevchenko const struct primes *p;
2222c64e9cbSAndy Shevchenko bool result;
2232c64e9cbSAndy Shevchenko
2242c64e9cbSAndy Shevchenko rcu_read_lock();
2252c64e9cbSAndy Shevchenko p = rcu_dereference(primes);
2262c64e9cbSAndy Shevchenko while (x >= p->sz) {
2272c64e9cbSAndy Shevchenko rcu_read_unlock();
2282c64e9cbSAndy Shevchenko
2292c64e9cbSAndy Shevchenko if (!expand_to_next_prime(x))
2302c64e9cbSAndy Shevchenko return slow_is_prime_number(x);
2312c64e9cbSAndy Shevchenko
2322c64e9cbSAndy Shevchenko rcu_read_lock();
2332c64e9cbSAndy Shevchenko p = rcu_dereference(primes);
2342c64e9cbSAndy Shevchenko }
2352c64e9cbSAndy Shevchenko result = test_bit(x, p->primes);
2362c64e9cbSAndy Shevchenko rcu_read_unlock();
2372c64e9cbSAndy Shevchenko
2382c64e9cbSAndy Shevchenko return result;
2392c64e9cbSAndy Shevchenko }
2402c64e9cbSAndy Shevchenko EXPORT_SYMBOL(is_prime_number);
2412c64e9cbSAndy Shevchenko
dump_primes(void)2422c64e9cbSAndy Shevchenko static void dump_primes(void)
2432c64e9cbSAndy Shevchenko {
2442c64e9cbSAndy Shevchenko const struct primes *p;
2452c64e9cbSAndy Shevchenko char *buf;
2462c64e9cbSAndy Shevchenko
2472c64e9cbSAndy Shevchenko buf = kmalloc(PAGE_SIZE, GFP_KERNEL);
2482c64e9cbSAndy Shevchenko
2492c64e9cbSAndy Shevchenko rcu_read_lock();
2502c64e9cbSAndy Shevchenko p = rcu_dereference(primes);
2512c64e9cbSAndy Shevchenko
2522c64e9cbSAndy Shevchenko if (buf)
2532c64e9cbSAndy Shevchenko bitmap_print_to_pagebuf(true, buf, p->primes, p->sz);
254*9ac17575SChristophe JAILLET pr_info("primes.{last=%lu, .sz=%lu, .primes[]=...x%lx} = %s\n",
2552c64e9cbSAndy Shevchenko p->last, p->sz, p->primes[BITS_TO_LONGS(p->sz) - 1], buf);
2562c64e9cbSAndy Shevchenko
2572c64e9cbSAndy Shevchenko rcu_read_unlock();
2582c64e9cbSAndy Shevchenko
2592c64e9cbSAndy Shevchenko kfree(buf);
2602c64e9cbSAndy Shevchenko }
2612c64e9cbSAndy Shevchenko
selftest(unsigned long max)2622c64e9cbSAndy Shevchenko static int selftest(unsigned long max)
2632c64e9cbSAndy Shevchenko {
2642c64e9cbSAndy Shevchenko unsigned long x, last;
2652c64e9cbSAndy Shevchenko
2662c64e9cbSAndy Shevchenko if (!max)
2672c64e9cbSAndy Shevchenko return 0;
2682c64e9cbSAndy Shevchenko
2692c64e9cbSAndy Shevchenko for (last = 0, x = 2; x < max; x++) {
2702c64e9cbSAndy Shevchenko bool slow = slow_is_prime_number(x);
2712c64e9cbSAndy Shevchenko bool fast = is_prime_number(x);
2722c64e9cbSAndy Shevchenko
2732c64e9cbSAndy Shevchenko if (slow != fast) {
274*9ac17575SChristophe JAILLET pr_err("inconsistent result for is-prime(%lu): slow=%s, fast=%s!\n",
2752c64e9cbSAndy Shevchenko x, slow ? "yes" : "no", fast ? "yes" : "no");
2762c64e9cbSAndy Shevchenko goto err;
2772c64e9cbSAndy Shevchenko }
2782c64e9cbSAndy Shevchenko
2792c64e9cbSAndy Shevchenko if (!slow)
2802c64e9cbSAndy Shevchenko continue;
2812c64e9cbSAndy Shevchenko
2822c64e9cbSAndy Shevchenko if (next_prime_number(last) != x) {
283*9ac17575SChristophe JAILLET pr_err("incorrect result for next-prime(%lu): expected %lu, got %lu\n",
2842c64e9cbSAndy Shevchenko last, x, next_prime_number(last));
2852c64e9cbSAndy Shevchenko goto err;
2862c64e9cbSAndy Shevchenko }
2872c64e9cbSAndy Shevchenko last = x;
2882c64e9cbSAndy Shevchenko }
2892c64e9cbSAndy Shevchenko
290*9ac17575SChristophe JAILLET pr_info("%s(%lu) passed, last prime was %lu\n", __func__, x, last);
2912c64e9cbSAndy Shevchenko return 0;
2922c64e9cbSAndy Shevchenko
2932c64e9cbSAndy Shevchenko err:
2942c64e9cbSAndy Shevchenko dump_primes();
2952c64e9cbSAndy Shevchenko return -EINVAL;
2962c64e9cbSAndy Shevchenko }
2972c64e9cbSAndy Shevchenko
primes_init(void)2982c64e9cbSAndy Shevchenko static int __init primes_init(void)
2992c64e9cbSAndy Shevchenko {
3002c64e9cbSAndy Shevchenko return selftest(selftest_max);
3012c64e9cbSAndy Shevchenko }
3022c64e9cbSAndy Shevchenko
primes_exit(void)3032c64e9cbSAndy Shevchenko static void __exit primes_exit(void)
3042c64e9cbSAndy Shevchenko {
3052c64e9cbSAndy Shevchenko free_primes();
3062c64e9cbSAndy Shevchenko }
3072c64e9cbSAndy Shevchenko
3082c64e9cbSAndy Shevchenko module_init(primes_init);
3092c64e9cbSAndy Shevchenko module_exit(primes_exit);
3102c64e9cbSAndy Shevchenko
3112c64e9cbSAndy Shevchenko module_param_named(selftest, selftest_max, ulong, 0400);
3122c64e9cbSAndy Shevchenko
3132c64e9cbSAndy Shevchenko MODULE_AUTHOR("Intel Corporation");
3142c64e9cbSAndy Shevchenko MODULE_LICENSE("GPL");
315